Stability Analysis of Chain, Mild and Passive Smoking Model

Abeer A. Alshareef, H. A. Batarfi^*

Department of Mathematics, Faculty of Science, King Abdulaziz University, Jeddah, Saudi Arabia.

DOI: 10.4236/ajcm.2020.101003                

Abstract

In this paper, the global stability of free smoking equilibrium point was evaluated and presented graphically. The linear stability of a developed mathematical model illustrates the effect on the population of chain, mild and passive smokers. MATLAB programming was used to simulate the solutions, the reproduction number R₀ and the nature of the equilibria.

Keywords

Smoking Model, Stability Analysis, Global Stability and Lyapunov Method, Qualitative Behavior, Passive Smokers

1. Introduction

Castillo et al. [1] proposed a simple mathematical model for giving up smoking. They assumed that the total population consists of three populations: potential smokers, i.e. people who do not smoke yet and may become smokers in the future (P), smokers (S) and people who have already quit smoking permanently (Q). Since then, many articles were published further developing this model. The modifications were produced either by adding a factor or a new population (see e.g. [2] [3] [4] [5]). However, in 2008, Sharomi and Gumel [6] modified the model by including a class of smokers who temporarily quit smoking. They studied the global stability of smoker-free equilibrium (SFE). In 2011, Zaman [7] derived the giving-up smoking model taking into account the occasional smoking population and presented its qualitative behavior. In 2014, Z. Alkhudhari et al. [8] modified the model by assuming that the smoking classes are divided into occasional and chain smokers and they studied the existence and stability of equilibrium points. In 2017, Matintu [9] proposed a new modification by considering the moderate, chain smokers, potential smokers, and temporarily quit smokers. In 2019, Ana et al. [10] constructed a new smoking model where the total population was partitioned into three classes: passive smoker $E_1$ (those at risk from others’ smoking), infected I (people who are addicted to tobacco and can pass on the habit) and $E_2$ (people who have stopped smoking but are at risk of relapse) described as:

${\stackrel{\dot{}}{E}}_1=\Delta -\frac{\beta }{N}{E}_{1}I-\left(\mu +v_1\right)E_1,$ (1a)

$\stackrel{\dot{}}{I}=\frac{\beta }{N}{E}_{1}I-\left(\mu +ϵ+\theta \right)I,$ (1b)

${\stackrel{\dot{}}{E}}_2=\theta I-\frac{\beta }{N}{E}_{2}I-\left(\mu +v_2+\omega \right)I,$ (1c)

where $E_1=E_1\left(t\right)$, $I=I\left(t\right)$ and $E_2=E_2\left(t\right)$. The parameters ( $\mu ,v_1,v_2,\omega ,\beta ,ϵ$ ) and $\theta$ represent respectively the rate of natural death, exit rate from $E_1$ to healthy population, exit rate form $E_2$ to healthy population, death rate of ex-smoker, infection rate of smoking, death rate of smokers and exit rate from I to $E_2$. $\Delta$ is the average number of healthy people who are at risk of becoming active smokers.

Our paper is divided into: section (2), evaluating the global stability of the model in [10]. The modification of the model was stated in section (3) to include two classes of smokers: chain smokers $S_1$ and mild smokers $S_2$, with the full study of stability analysis and effect of each group on the population behavior (qualitative behavior).

2. The Global Stability of Equilibria

The equilibria of smoking-free equilibrium point in system (1a)-(1b), are given in [10] as

$P_1=\left(E_{1,0},I_0,E_{2,0}\right)=\left(\frac{\Delta }{\mu +v_2},0,0\right),$ (2)

and the basic reproduction number ${\mathcal{R}}_0$ is defined as ${\mathcal{R}}_0=\frac{\beta }{\mu +ϵ+\theta }$.

Theorem 1.

If $\beta \le {\delta }_4$ and ${\mathcal{R}}_0\le 1$, then $P_1$ is globally asymptotically stable.

Proof.

The Lyapunov function is defined as

$V_1=E_{1,0}\left(\frac{E_1}{E_{1,0}}-1-\ln \frac{E_1}{E_{1,0}}\right)+E_2+I,$

hence,

${\stackrel{\dot{}}{V}}_1=\left(1-\frac{E_{1,0}}{E_1}\right)+\stackrel{\dot{}}{I}+{\stackrel{\dot{}}{E}}_2,$ (3)

By substituting (1a)-(1b) in (3),

$\begin{array}{c}{\stackrel{\dot{}}{V}}_1=\left(1-\frac{E_{1,0}}{E_1}\right)\left(\Delta -\frac{\beta }{N}{E}_{1}I-\left(\mu +v_1\right)E_1\right)+\frac{\beta }{N}I\left(E_1+E_2\right)\\ -\left(\mu +ϵ+\theta \right)I+\theta I-\frac{\beta }{N}{E}_{2}I-\left(\mu +v_2+\omega \right)E_2,\end{array}$

$\begin{array}{c}{\stackrel{\dot{}}{V}}_1=\left(1-\frac{E_{1,0}}{E_1}\right)\left(\Delta -\left(\mu +v_1\right)E_1\right)-\frac{\beta }{N}{E}_{1}I+\frac{\beta }{N}{E}_{0,1}I+\frac{\beta }{N}I\left(E_1+E_2\right)\\ -\left(\mu +ϵ+\theta \right)I+\theta I-\frac{\beta }{N}{E}_{2}I-\left(\mu +v_2+\omega \right)E_2,\end{array}$

Since $E_{1,0}=N$, then we have

${\stackrel{\dot{}}{V}}_1=\left(1-\frac{E_{1,0}}{E_1}\right)\left(\Delta -{\delta }_1{E}_1\right)-{\delta }_2{E}_2+\left(\beta -\left(\mu +ϵ\right)\right)I.$

where ${\delta }_1=\mu +v_1$ and ${\delta }_2=\mu +v_2+\omega$. Using (2) we get

${\stackrel{\dot{}}{V}}_1=-{\delta }_1\frac{{\left(E_{1,0}-E_1\right)}^2}{E_1}-{\delta }_2{E}_2+\left(\beta -\left(\mu +ϵ\right)\right)I.$

By re-arranging the equation, we get

${\stackrel{\dot{}}{V}}_1=-{\delta }_1\frac{{\left(E_{1,0}-E_1\right)}^2}{E_1}-{\delta }_2{\mathcal{R}}_0{E}_2+\frac{{\delta }_2}{{\delta }_3}\left({\delta }_3\left({\mathcal{R}}_0-1\right)\right)E_2+\left(\beta -{\delta }_4\right)I,$ (4)

where ${\delta }_3=\mu +ϵ+\theta$ and ${\delta }_4=\mu +ϵ$.

It is clear that ${\stackrel{\dot{}}{V}}_1\le 0$ when $\beta \le {\delta }_4$ and ${\mathcal{R}}_0\le 1$ for all $E_1,I,E_2>0$. Hence, the solutions of systems (1) is limited to $\Omega$, the largest invariant subset of ${\stackrel{\dot{}}{V}}_1=0$, where $\Omega$ is the region of solutions given in [10]. From (4), we see that ${\stackrel{\dot{}}{V}}_1=0$ if and only if $E_1=E_{1,0}$ and $I=E_2=0$. By using LaSalle’s invariance principle, the point $P_1$ is globally asymptotically stable.

□

3. The Modified Model of Smoking

In this section, system (1) was modified to include two classes of smokers, chain smokers $S_1$ and mild smoker $S_2$ (see Figure 1). Accordingly, our model of ordinary differential equations (ODEs) are

Figure 1. Graphic presentation of system (6).

${\stackrel{\dot{}}{E}}_1\left(t\right)=\text{Δ}-\frac{\beta }{N}{E}_1\left(t\right)\left(\eta S_1\left(t\right)+S_2\left(t\right)\right)-\left(\mu +v_1\right)E_1\left(t\right),$ (5a)

${\stackrel{\dot{}}{S}}_1\left(t\right)=a\frac{\beta }{N}\left(\eta S_1\left(t\right)+S_2\left(t\right)\right)\left(E_1\left(t\right)+E_2\left(t\right)\right)-\left(\mu +{ϵ}_1+{\theta }_1\right)S_1\left(t\right),$ (5b)

${\stackrel{\dot{}}{S}}_2\left(t\right)=\left(1-a\right)\frac{\beta }{N}\left(\eta S_1\left(t\right)+S_2\left(t\right)\right)\left(E_1\left(t\right)+E_2\left(t\right)\right)-\left(\mu +{ϵ}_2+{\theta }_2\right)S_2\left(t\right),$ (5c)

${\stackrel{\dot{}}{E}}_2\left(t\right)={\theta }_1{S}_1\left(t\right)+{\theta }_2{S}_2\left(t\right)-\frac{\beta }{N}{E}_2\left(t\right)\left(\eta S_1\left(t\right)+S_2\left(t\right)\right)-\left(\mu +v_2+\omega \right)E_2\left(t\right),$ (5d)

where the parameters $\beta ,\omega ,v_1,v_2$ and $\mu$ are defined as in the system (1a)-(1c), ${ϵ}_1$ and ${ϵ}_2$ are death rate of chain smokers $S_1$ and mild smokers $S_2$ respectively, while parameters ${\theta }_1$ and ${\theta }_2$ are exit rates from chain smokers and mild smokers to the healthy population (outside population N). We assume that the exposed people become either a chain or mild smoker at probabilities (1 − a) and a with $0<a<1$ ). The chain smokers have a higher probability of generating more new smokers (by a factor $\eta \ge 1$ relative to the mild smokers). A simplified form of the model (5) can be written as

${\stackrel{\dot{}}{E}}_1\left(t\right)=\text{Δ}-\frac{\beta }{N}{E}_1\left(t\right)\left(\eta S_1\left(t\right)+S_2\left(t\right)\right)-d_1{E}_1\left(t\right),$ (6a)

${\stackrel{\dot{}}{S}}_1\left(t\right)=a\frac{\beta }{N}\left(\eta S_1\left(t\right)+S_2\left(t\right)\right)\left(E_1\left(t\right)+E_2\left(t\right)\right)-d_2{S}_1\left(t\right),$ (6b)

${\stackrel{\dot{}}{S}}_2\left(t\right)=\left(1-a\right)\frac{\beta }{N}\left(\eta S_1\left(t\right)+S_2\left(t\right)\right)\left(E_1\left(t\right)+E_2\left(t\right)\right)-d_3{S}_2\left(t\right),$ (6c)

${\stackrel{\dot{}}{E}}_2\left(t\right)={\theta }_1{S}_1\left(t\right)+{\theta }_2{S}_2\left(t\right)-\frac{\beta }{N}{E}_2\left(t\right)\left(\eta S_1\left(t\right)+S_2\left(t\right)\right)-d_4{E}_2\left(t\right),$ (6d)

where $d_1=\mu +v_1$, $d_2=\left(\mu +{ϵ}_1+{\theta }_1\right)$, $d_3=\left(\mu +{ϵ}_2+{\theta }_2\right)$ and $d_4=\left(\mu +v_2+\omega \right)$. The initial conditions of system (6), are given by

$E_1\left(\alpha \right)={\phi }_1\left(\alpha \right)$, $S_1\left(\alpha \right)={\phi }_2\left(\alpha \right)$, $S_2\left(\alpha \right)={\phi }_3\left(\alpha \right)$, $E_2\left(\alpha \right)={\phi }_4(\alpha )$

where $\alpha \in R$, $\left({\phi }_1\left(\alpha \right),{\phi }_2\left(\alpha \right),{\phi }_3\left(\alpha \right),{\phi }_4\left(\alpha \right)\right)\in C{R}_{+}^4$ and C is the Banach space of continuous functions mapping the interval $\left(-\infty ,0\right]$ into $R_{+}^4$. By the fundamental theory of functional differential equations [11], system (6) has a unique solution satisfying the initial conditions above.

3.1. Non-Negativity and Boundedness of Solutions

Theorem 2. For system (6), there exist a positive number K such that the compact set $\text{Λ}$,

$\text{Λ}=\left\{\left(E_1,S_1,S_2,E_2\right)\in R_{+}^4:0\le E_1,S_1,S_2,E_2\le K\right\}.$

is positively invariant.

Proof. We have

${{\stackrel{\dot{}}{E}}_1\left(t\right)|}_{E_1=0}=\text{Δ}>0$,

${{\stackrel{\dot{}}{S}}_1\left(t\right)|}_{S_1=0}=a\frac{\beta }{N}{S}_2\left(E_1+E_2\right)>0$ for   $S_2,E_1$ and $E_2>0$,

${{\stackrel{\dot{}}{S}}_2\left(t\right)|}_{S_2=0}=\left(1-a\right)\frac{\beta }{N}\eta S_1\left(E_1+E_2\right)>0$ for   $S_1,E_1$ and $E_2>0$,

${{\stackrel{\dot{}}{E}}_2\left(t\right)|}_{E_2=0}={\theta }_1{S}_1+{\theta }_2{S}_2>0$ for   $S_1$ and   $S_2>0$,

so the solutions are non-negative. Now, to prove the boundedness of the solutions, we let $N\left(t\right)=E_1\left(t\right)+S_1\left(t\right)+S_2\left(t\right)+E_2\left(t\right)$, hence

$\begin{array}{c}\stackrel{\dot{}}{N}=\text{Δ}-\frac{\beta }{N}{E}_1\left(\eta S_1+S_2\right)-d_1{E}_1+a\frac{\beta }{N}\left(\eta S_1+S_2\right)\left(E_1+E_2\right)-d_2{S}_1-d_3{S}_2\\ +\left(1-a\right)\frac{\beta }{N}\left(\eta S_1+S_2\right)\left(E_1+E_2\right)+{\theta }_1{S}_1+{\theta }_2{S}_2-\frac{\beta }{N}{E}_2\left(\eta S_1{S}_2\right)-d_4{E}_2\\ \le \text{Δ}-d_1{E}_1-d_5{S}_1-d_6{S}_2-d_4{E}_2\le \text{Δ}-\bar{d}N\left(t\right),\end{array}$

where $d_5=\mu +{ϵ}_1$, $d_6=\mu +{ϵ}_2$ and $\bar{d}=\min \left\{d_i\right\}$, $i=1,4,5,6$. That is, if $N\left(0\right)\le 0$ then $N\left(t\right)<K$ where $K=\text{Δ}/\bar{d}$. It follows that if $E_1\left(0\right)+S_1\left(0\right)+S_2\left(0\right)+E_2\left(0\right)\le K$ then $E_1\left(t\right),S_1\left(t\right),S_2\left(t\right),E_2\left(t\right)\le K$.

□

3.2. The Equilibria and Basic Reproduction Number ${\mathcal{R}}_0$

Theorem 3. 1) If ${\mathcal{R}}_0\le 1$, then there exists only one equilibrium point called smoking-free equilibrium $P_0$.

2) If ${\mathcal{R}}_0>1$, there exist two equilibria that are smoking-free equilibrium $P_0$ and smoking present equilibrium $P_1$.

Proof.

The system (6) has two equilibrium points which are: Smoking-free equilibrium point

$P_0=\left(E_{1,0},S_{1,0},S_{2,0},E_{2,0}\right)=\left(\frac{\text{Δ}}{d_1},0,0,0\right),$ (7)

and smoking present equilibrium point

$P_1=\left(E_1^{\text{*}},S_1^{\text{*}},S_2^{\text{*}},E_2^{\text{*}}\right)$

where

$E_1^{*}=\frac{{\mathcal{R}}_0\text{Δ}}{\left({\mathcal{R}}_0-1\right)\beta +{\mathcal{R}}_0\left(\mu +v_1\right)},$

$S_1^{*}=a\frac{\beta \left({\mathcal{R}}_0-1\right)}{{\mathcal{R}}_0{d}_2}\left(x_1+E_2^{*}\right),$

$S_2^{*}=\left(1-a\right)\frac{\beta \left({\mathcal{R}}_0-1\right)}{{\mathcal{R}}_0{d}_3}\left(x_1+E_2^{*}\right),$

$E_2^{*}=\frac{\text{Δ}\beta \left({\mathcal{R}}_0-1\right)\left[\left(1-a\right){\theta }_2{d}_2+a{\theta }_1{d}_3\right]}{d_2{d}_3{x}_2\left[\beta \left({\mathcal{R}}_0-1\right)+{\mathcal{R}}_0{d}_1\right]},$

and

$x_1=\left(\frac{\text{Δ}{\mathcal{R}}_0}{{\mathcal{R}}_0{d}_1+\beta \left({\mathcal{R}}_0-1\right)}\right),$

$x_2=d_4+\left(\frac{\beta \left({\mathcal{R}}_0-1\right)}{{\mathcal{R}}_0{d}_2{d}_3}\right)\left(d_2{d}_6\left(1-a\right)+a{d}_3{d}_5\right).$

Using the next generation method [12], the matrices F and V are

$F=\left(\begin{array}{cccc}0& -\eta \beta & -\beta & 0\\ 0& a\eta \beta & a\beta & 0\\ 0& \left(1-a\right)\eta \beta & \left(1-a\right)\beta & 0\\ 0& 0& 0& 0\end{array}\right),V=\left(\begin{array}{cccc}{d}_1& 0& 0& 0\\ 0& d_2& 0& 0\\ 0& 0& d_3& 0\\ 0& 0& 0& d_4\end{array}\right),$

and hence evaluate the reproduction number ${\mathcal{R}}_0$ which is the largest eigenvalue of $F{\left(V\right)}^{-1}$,

${\mathcal{R}}_0=\beta \frac{\left(1-a\right)d_2+a\eta d_3}{d_2{d}_3}.$ (8)

□

3.3. Local Stability of Smoker Free Equilibrium

Theorem 4.

If ${\mathcal{R}}_0<1$, then $P_0$ is locally asymptotically stable.

Proof. The Jacobian matrix of system (6) at $P_0$ is

$J\left(P_0\right)=\left(\begin{array}{cccc}-d_1& -\eta \beta & \beta & 0\\ 0& \eta a\beta -d_2& a\beta & 0\\ 0& \left(1-a\right)\beta \eta & \left(1-a\right)\beta -d_3& 0\\ 0& {\theta }_1& {\theta }_2& -d_4\end{array}\right),$

with the eigenvalues ${\lambda }_1=-d_1$, ${\lambda }_2=-d_4$, and ${\lambda }_{3,4}=c\pm b$, where

$c=\left(\frac{1}{2}\left[\beta \left(1-a\right)-d_3\right]+\left[\beta a\eta -d_2\right]\right),$

and

$b=\frac{1}{2}\sqrt{2\left({\mathcal{R}}_0-1\right)+{\left(d_2-a\beta \eta \right)}^2+{\left(d_3-\left(1-a\right)\beta \right)}^2+2a\eta \beta \left(\beta \left(1-a\right)\right)}.$

It is easy to show that $c<0$ when ${\mathcal{R}}_0<1$. To study the nature of b, we start with the fact that,

$\beta \left(\left(1-a\right)d_2+a\eta d_3\right)<d_2{d}_3,$

hence,

$\beta \left(1-a\right)d_2<d_2{d}_3$, and $\beta a\eta d_3<d_2{d}_3$,

$\beta \left(1-a\right)<d_3$, and $\beta a\eta <d_2$.

Next, we use the numerical approach by assuming two small positive real numbers ${\xi }_1$ and ${\xi }_2$ such that $\beta \eta a+{\xi }_1=d_2$ and $\beta \left(1-a\right)+{\xi }_2=d_3$ hence

$\begin{array}{c}b=\frac{1}{2}\sqrt{2\left({\mathcal{R}}_0-1\right)+{\xi }_1^2+{\xi }_2^2+2\left(d_2-{\xi }_1\right)\left(d_3-{\xi }_2\right)}\\ =\frac{\sqrt{2}}{2}\sqrt{\left({\mathcal{R}}_0-1\right)+\frac{1}{2}{\xi }_1^2+\frac{1}{2}{\xi }_2^2+d_2{d}_3-{\xi }_2{d}_2-{\xi }_1{d}_3+{\xi }_1{\xi }_2}\\ =\frac{\sqrt{2}}{2}\sqrt{\left({\mathcal{R}}_0-1\right)+\frac{1}{2}{\left({\xi }_1+{\xi }_2\right)}^2+d_2{d}_3-{\xi }_2{d}_2-{\xi }_1{d}_3}.\end{array}$

Since $\frac{1}{2}{\left({\xi }_1+{\xi }_2\right)}^2\cong 0$ and $d_2{d}_3<{\xi }_1{d}_2+{\xi }_2{d}_3$ then $d_2{d}_3<d_2+d_3$ and

$\left({\mathcal{R}}_0-1\right)+\frac{1}{2}{\left({\xi }_1+{\xi }_2\right)}^2+d_2{d}_3-{\xi }_2{d}_2-{\xi }_1{d}_3<0.$

Thus ${\lambda }_{3,4}$ are complex numbers. Finally, we see that all the real parts of ${\lambda }_i,i=1,2,3,4$ are negative so the point $P_0$ is locally asymptotically stable.

□

3.4. The Global Stability of Smoker Free Equilibrium

Theorem 5. If $\beta \eta \le d_5,\beta \le d_6$ and ${\mathcal{R}}_0\le 1$, then $P_0$ is globally asymptotically stable.

Proof. As done in Theorem (1), we defined the Lyapunov function as

$V_2=E_{1,0}\left(\frac{E_1}{E_{1,0}}-1-\ln \frac{E_1}{E_{1,0}}\right)+S_1+S_2+E_2,$

hence

${\stackrel{\dot{}}{V}}_2=\left(1-\frac{E_{1,0}}{E_1}\right){\stackrel{\dot{}}{E}}_1+{\stackrel{\dot{}}{S}}_1+{\stackrel{\dot{}}{S}}_2+{\stackrel{\dot{}}{E}}_2.$ (9)

By substituting (6) in (9),

$\begin{array}{c}{\stackrel{\dot{}}{V}}_2=\left(1-\frac{E_{1,0}}{E_1}\right)\left(\text{Δ}-\frac{\beta }{N}{E}_1\left(\eta S_1+S_2\right)-d_1{E}_1\right)+\frac{a\beta }{N}\left(\eta S_1+S_2\right)\left(E_1+E_2\right)-d_2{S}_1\\ +\frac{\left(1-a\right)\beta }{N}\left(\eta S_1+S_2\right)\left(E_1+E_2\right)-d_3{S}_2+{\theta }_1{S}_1+{\theta }_2{S}_2-\frac{\beta }{N}{E}_2\left(\eta S_1+S_2\right)-d_4{E}_2\\ =\left(1-\frac{E_{1,0}}{E_1}\right)\left(\text{Δ}-d_1{E}_1\right)+\frac{\beta }{N}{E}_{1,0}\left(\eta S_1+S_2\right)-d_2{S}_1-d_3{S}_2+{\theta }_1{S}_1+{\theta }_2{S}_2-d_4{E}_2\end{array}$

By using the equilibrium conditions $\text{Δ}=d_1{E}_{1,0}$ and $E_{1,0}=N$,

$\begin{array}{c}{\stackrel{\dot{}}{V}}_2=-d_1\left(\frac{{\left(E_1-E_{1,0}\right)}^2}{E_1}\right)+\left(\beta \eta -d_5\right)S_1+\left(\beta -d_6\right)S_2-d_4{E}_2\\ =-d_1\left(\frac{{\left(E_1-E_{1,0}\right)}^2}{E_1}\right)+\left(\beta \eta -d_5\right)S_1+\left(\beta -d_6\right)S_2\\ +\frac{d_4}{c_1}\left(c_1\left[{\mathcal{R}}_0-1\right]\right)E_2-{\mathcal{R}}_0{d}_4{E}_2,\end{array}$ (10)

where $c_1=d_2{d}_3$ and $c_2=\beta \left(1-a\right)d_2+a\eta d_3$. Clearly ${\stackrel{\dot{}}{V}}_2\le 0$, if $\beta \eta \le d_5$, $\beta \le d_6$ and ${\mathcal{R}}_0\le 1$, for all $E_1,S_1,S_2,E_2>0$. Hence, the solutions of system (6) are limited to $\omega$, the largest invariant subset of ${\stackrel{\dot{}}{V}}_2=0$. From (10) we see that ${\stackrel{\dot{}}{V}}_2=0$ if and only if $E_1=E_{1,0}$ and $S_1=S_2=E_2=0$. Using LaSalle’s invariance principle, we show that $P_0$ is globally asymptotically stable.

□

3.5. Numerical Simulations and Discussion

Using MATLAB programming, the solution of system (6) is evaluated by assuming values for the parameters as $\text{Δ}=50$, ${ϵ}_1=0.03$, ${ϵ}_2=0.02$, ${\theta }_1=0.4$, ${\theta }_2=0.5$, $\mu =0.05$, $v_1=0.2$, $v_2=0.3$, $N=200$, and $\omega =0.001$ with two sets of values for $\beta ,\eta$ and a.

Set (1): ${\mathcal{R}}_0\le 1$, $\beta =0.05$, $\eta =1.5$, and $a=0.4$.

Set (2): ${\mathcal{R}}_0>1$, $\beta =0.5$, $\eta =2.5$, and $a=0.6$.

A three sets of values are applied and tested for the initial conditions:

IC1: $E_1\left(0\right)=890,S_1\left(0\right)=120,S_2\left(0\right)=90,E_2\left(0\right)=40$.

IC2: $E_1\left(0\right)=480,S_1\left(0\right)=60,S_2\left(0\right)=40,E_2\left(0\right)=30$.

IC3: $E_1\left(0\right)=290,S_1\left(0\right)=30,S_2\left(0\right)=20,E_2\left(0\right)=10$.

In Figures 2-5, the behavior of each smoking population is plotted as time increases. Figure 2 shows the fast decrease of passive smokers for both ${\mathcal{R}}_0\le 1$ and ${\mathcal{R}}_0>1$ within a short period of time and reaches or approaches the steady

Figure 2. The plot of the passive smokers ( $E_1$ ), for set (1), ${\mathcal{R}}_0<1$ and set (2), ${\mathcal{R}}_0>1$ together with IC1 (dash), IC2 (solid) and IC3 (dash dot).

Figure 3. The plot of the chain smokers ( $S_1$ ), for set (1), ${\mathcal{R}}_0<1$ and set (2), ${\mathcal{R}}_0>1$ together with IC1 (dash), IC2 (solid) and IC3 (dash dot).

Figure 4. The plot of the mild smokers ( $S_2$ ), for set (1), ${\mathcal{R}}_0<1$ and set (2), ${\mathcal{R}}_0>1$ together with IC1 (dash), IC2 (solid) and IC3 (dash dot).

Figure 5. The plot of the ex-smokers ( $E_2$ ), for set (1), ${\mathcal{R}}_0<1$ and set (2), ${\mathcal{R}}_0>1$ together with IC1 (dash), IC2 (solid) and IC3 (dash dot).

state. Similar behavior was observed on the chain and mild populations for ${\mathcal{R}}_0\le 1$, Figure 3 and Figure 4 which reached zero as time increased. For ${\mathcal{R}}_0>1$, a peak appears at the beginning of the time and decreases to the steady state as time increases. This behavior was seen on the ex-smokers’ population for both ${\mathcal{R}}_0\le 1$ and ${\mathcal{R}}_0>1$. The choice of initial values of the three sets has no effect on the behaviour of the different population types. It tends to zero for the chain, mild and ex-smokers for ${\mathcal{R}}_0<1$ as time increases. $E_2$ behaves similarly for any initial conditions, where we observe that $E_2$ tends to zero ( ${\mathcal{R}}_0\le 1$ ) or steady state value ( ${\mathcal{R}}_0>1$ ). Similarly for other types of population which means that $P_0$ and $P_1$ are globally asymptomatically stable. Figure 6 shows the phase plane space plots that explain the relationships between the different populations. These relationships are the same for different initial conditions. Next, we consider the relationships of ex-smokers’ population versus the chain/mild smokers

Figure 6. The phase plane plot of the populations of system (6) using the values of the parameters of set (1) with the initial conditions $E_1\left(0\right)=30$, $S_1\left(0\right)=290$, $S_2\left(0\right)=20$ and $E_2\left(0\right)=10$.

$S_1$ and $S_2$. As the chain smokers’ number decreases, we notice an increase on the population of ex-smokers ( $E_2$ ) up to a maximum value at around one third of chain smokers’ initial values and around two-thirds of the mild smokers’ initial value, Figure 5(c) and Figure 5(d). Finally, Figure 5(e) illustrates the impact of decreasing the number of chain smokers which leads to a decrease in the number of mild smokers.

As we observe from the Figures, we stress the importance of decreasing the number of chain smokers (either by using the media or educational direct seminars) on increasing the number of passive smokers (as the ex-smokers and mild are both joined) to create a healthy action.

Finally, we observe that when the contact between the two classes of smokers with the other populations is weak (small β) then the effect of the smoker’s classes (chain/mild) was not being enough to generate the new individuals to the smokers’ classes. So, the number of smoker’s groups decreased and approached to zero (because the death or stop smoking). On the other hand, if β increase, we found that the number of smoker’s groups increased also, since the rate of interaction between them and other groups (or the effect by smokers) was growing. However, we can see that, one of the solutions which is introduced to reduce the spread of smoking is prevent it in the larger places (such as educations) in order to minimize the rate of contact to be very small or zero.

4. Conclusion

The modified system/model shows interesting behavior between the different types of smoker’s populations. The main advantage was the more decrease in the size of heavy/mild smokers in comparison with earlier studies. This should encourage others (ex- or passive) smokers to interact better with chain/mild smokers to increase the number of quitters. By comparing our results with the results in [10], we see that they are approximately similar to giving the behavior of each group, and the difference is that, our model explains and shows the effect of both of the chain smokers and mild smokers separately on the exposed population. Hence, we get that the chain smokers affect more than the mild smokers on the exposed groups. So, our model is more real.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1]	Castillo, C., Jordan, G. and Rodriguez, A. (1997) Mathematical Models for the Dynamics of Tobacco Use, Recovery, and Relapse. Technical Report Series.
[2]	Awan, A., Sharif, A., Hussain, T. and Ozair, M. (2017) Smoking Model with Cravings to Smoke. Advanced Studies in Biology, 9, 31-41. https://doi.org/10.12988/asb.2017.61245
[3]	Alkhudhari, Z., Al-Sheikh, S. and Al-Tuwairqi, S. (2014) Global Dynamics of a Mathematical Model on Smoking. Applied Mathematics, 2014, Article ID: 847075. https://doi.org/10.1155/2014/847075
[4]	Sallew, G. (2016) Mathematical Analysis on the Spread and Control of Global Dynamics of Tobacco Smoking with Induced Death. Journal of Mathematics and Computer, 3, 15-26.
[5]	Lahrouz, A., Omari, L., Kiouach, D. and Belmaati, A. (2011) Deterministic and Stochastic Stability of a Mathematical Model of Smoking. Statistics and Probability Letters, 81, 1276-1284. https://doi.org/10.1016/j.spl.2011.03.029
[6]	Sharomi, O. and Gumel, A. (2008) Curtailing Smoking Dynamics: A Mathematical Modeling Approach. Applied Mathematics and Computation, 195, 475-499. https://doi.org/10.1016/j.amc.2007.05.012
[7]	Zaman, G. (2011) Qualitative Behavior of Giving up Smoking Models. Bulletin of the Malaysian Mathematical Sciences Society, 34, 403-415.
[8]	Alkhudhari, Z., Al-Sheikh, S. and Al-Tuwairqi, S. (2014) The Effect of Occasional Smokers on the Dynamics of a Smoking Model. International Mathematical Forum, 9, 1207-1222. https://doi.org/10.12988/imf.2014.46120
[9]	Matintu, S.A. (2017) Smoking as Epidemic: Modeling and Simulation Study. American Journal of Applied Mathematics, 5, 31-38. https://doi.org/10.11648/j.ajam.20170501.14
[10]	Pulecio-Montoya, A., Lopez-Monteneqro, L. and Benavides, L. (2019) Analysis of a Mathematical Model of Smoking. Contemporary Engineering Sciences, 12, 117-129. https://doi.org/10.12988/ces.2019.9517
[11]	Hale, J. and Lunel, S. (1993) Introduction to Functional Differential Equations. Science and Business Media. https://doi.org/10.1007/978-1-4612-4342-7
[12]	Heffernan, J., Smith, R. and Wahl, L. (2005) Perspectives on the Basic Reproduction Ratio. Journal of the Royal Society Interface, 2, 281-293. https://doi.org/10.1098/rsif.2005.0042