A Deterministic Mathematical Dynamic System Based on the PSITPS Model for Modeling the COVID-19 Epidemic

In 1760, Daniel Bernoulli proposed a resolution to his mathematical model pertaining to smallpox [1], and modeling infectious diseases are still as important today. The use of mathematical models to research epidemic causes and effects, however, did not gain widespread acceptance until the 20th century. The first model was created in 1927 by Kermack and McKendrick and is known as the SIR model [2]. This model classifies the population into three distinct groups: susceptible, infected, and removed. The model takes into account group interactions that are based on rates of disease transmission and eradication as well as a constant total population.

The individuals that have been removed can be employed to transform a population size that is dependent on time into a population size that is constant, even though in some models they are used to represent vaccinated individuals. There are advanced modifications of the SIR model, such as the SIRD, SEIR, and SEIRD models [3-5]. Hundreds of viruses comprise the coronavirus family. The SARS-CoV, or severe acute respiratory syndrome coronavirus, was identified in November 2002. In September 2012, researchers identified the Middle East respiratory syndrome coronavirus (MERS-CoV). Both were linked to serious respiratory conditions and high mortality rates [6, 7]. The identification of SARS-CoV-2, a novel coronavirus associated with severe acute respiratory syndrome, has led to the emergence of a contagious disease often referred to as coronavirus disease 2019 (COVID-19) [8].

In Wuhan, China, the virus was initially identified and reported in December 2019. SARS-CoV-2 is a very contagious virus and rapidly spreads around the world [9]. More than 2.1 million COVID-19 cases have been confirmed by the WHO, together with 142,229 fatalities [10]. Iraq is one of the nations that has had significant repercussions as a result of the viral outbreak. On February 25, 2020, the first recorded case of COVID-19 was confirmed in the Najaf governorate for an Iranian student who came from Iran [11], supported by four cases from one family in the Kirkuk governorate on February 26, 2020. All of these cases had a travel history to Iran [12].

COVID-19 cases in total from January 3, 2020, to December 8, 2022, is 2,463,296 confirmed cases with 25,363 deaths in Iraq [13]. While in Najaf Governorate, the cumulative total of COVID-19 cases from March 31, 2020, to May 9, 2022, is 95,147 confirmed cases with 728 deaths [14]. However, according to available data, the cumulative count of COVID-19 cases in Najaf from January 1, 2021, to December 31, 2021, amounted to 68,394, resulting in 397 fatalities [15].

Many studies have been performed on the dynamics of COVID-19.

Abidemi et al. [16] offered a carefully designed deterministic compartmental model to explore the influence of several pharmacological and non-pharmacological control methods on COVID-19 population dynamics in Malaysia. The influence of different control combination techniques incorporating the use of safeguarding oneself, tracking of contacts and laboratory tests, and treatment control measures on disease transmission is demonstrated using numerical simulations. According to numerical simulations, the use of each technique analysed can dramatically lower both the prevalence and incidence of COVID-19 within the general community.

Lobato et al. [17] proposed a technique for parameter estimation in compartmental models through the use of a dynamic data segmentation approach. The likelihood of transmission by touch for each wave is represented by a time-dependent function. The findings offer details on the epidemic's features for each nation as well as input for economic and governmental decision-making.

Rida et al. [18] presented a Mathematical SEIAR Model to analyse the transmission dynamics of the ongoing COVID-19 epidemic. Where the population was categorized into five distinct groups, namely susceptible $S(t)$, exposed $E(t)$, symptomatic $I(t)$, asymptomatic $A(t)$, and recovered $R(t)$. The reproduction number $\mathcal{R}_0$ is used to conduct a sensitivity analysis and identify the most important variables from which to calculate the COVID-19 propagation volume. A numerical simulation was given and is based on data to estimate COVID19's worldwide spread.

Mustafa and Ahmad [19] investigated the impact of contaminated objects on a SIRS model, including compartments for vaccination and hospitalization. The model's solutions are shown to be positive and bounded. When the value of the Basic Reproduction Number is below 1, the examination of the Disease-Free Equilibrium point demonstrates its worldwide stability. Studies on the presence and distinctiveness of the Endemic Equilibrium point are conducted when the Basic Reproduction Number surpasses one.

Hamou et al. [20] studied the mathematical modeling of undiscovered new coronavirus infections in Morocco. The model's parameters and undiscovered COVID-19 instances were evaluated using the least squares approach and were discarded. The researchers employed the MATLAB programming language to explore unexplored scenarios in Morocco and confirm projected results.

Salam et al. [21] developed the SIUWR model they previously used and added certain transmission parameters. The model has five compartments called susceptible people S, people who are infected but have no symptoms I, unreported cases of symptomatic infection U, reported infected people with symptoms W, and recovered people R. The calculation of the Basic Reproduction Number and elastic characteristic coefficients was performed at the Equilibrium Points. The stated model's validity was tested by comparing it to some real data that was collected in Iraq and France between January 1^st and December 25^th, 2021. The results offered some biological interpretations, which may be utilized to help control this epidemic.

Rai [22] proposed a compartmental SVEIDRS Epidemic Model in order to examine the diffusion dynamics associated with COVID-19. The whole population (N) has been split into six compartments in the bio-mathematical model used. S(t) Amount of people who have not been exposed to the virus and have not developed immunity to it. V(t) is the number of people who were immunized against the virus. E(t) Amount of people who have been exposed to the virus but have not exhibited symptoms and are still in the incubation stage. I(t) Amount of people that have previously been exposed to the virus, are exhibiting illness symptoms and can spread the sickness. D(t) Amount of people who died as a result of complications induced by the virus. R(t) Number of people who have recovered from the virus and are immune to it for a set period of time. Examined are the results of vaccination as well as its effectiveness. A sensitivity analysis is conducted to examine the impact of model parameters on the propagation of a disease in the presence of vaccination. The possibility of COVID-19 is decreased, and the spread of the disease is controlled through the use of efficient vaccinations, improvements to healthcare infrastructure, and intentional measures to minimize human-to-human contact.

In a recent study, Mohammed and Khoshnaw [23] developed a mathematical model that incorporates a vaccination compartment. The model employed in this study adopts an equation-based methodology and aims to offer a full comprehension of the dynamics pertaining to vaccination. The present investigation centers on the analysis of empirical evidence derived from officially confirmed cases in the Kurdistan Region of Iraq during the period spanning from July 17, 2021, to January 1, 2022. In this study, the authors utilized MATLAB's computational tools to compare the total number of infected individuals with both model predictions and empirical data. The computational analysis reveals a striking resemblance between the dynamics of model outcomes and the actual, validated instances. The examination of elasticity coefficients provides a comprehensive comprehension of the impact of vaccination on transmission. The results of this model possess the capacity to provide a valuable contribution to the worldwide effort of effectively managing this ailment by offering supplementary recommendations and enhancements.

The authors, Satar and Naji [24] conducted a study whereby they introduced and examined a mathematical framework that effectively represents the transmission dynamics of COVID-19 throughout the population of asymptomatic individuals. A thorough examination was undertaken to analyse every component of the system's solution. A computation was conducted to determine the Equilibrium Points and the Basic Reproduction Number. The researchers performed local stability studies for both the equilibrium state without disease and the equilibrium state with disease present. In contrast to the Castillo-Chavez theorem, which was utilised to examine the global stability of the Disease-Free Equilibrium, a geometric approach was adopted to analyse the global dynamics of the Endemic Equilibrium. The analysis uncovered that the system's transcritical bifurcation occurs when there is an absence of sickness. The confirmation of the findings and the determination of the mechanisms governing the transmission of the disease were accomplished through the use of a numerical simulation, specifically applied to the population of Iraq.

To our knowledge, no model examines the effect of protection (vaccination) after treatment. According to what was mentioned above, we extended the model introduced by Wameko et al. [25] to the PSITPS Model for the COVID-19 epidemic and demonstrated the effect of protection (vaccination) after treatment in Najaf.

The PSITPS Model, which emphasizes post-treatment protection through vaccination, was examined in the context of Najaf, Iraq. Real data from the year 2021, obtained from the Public Health Department and the Statistics Division of the Training and Human Development Centre of the Najaf Health Directorate, Ministry of Health, were utilized to conduct numerical simulations and obtain the corresponding results.

This paper is arranged as follows, in the second Section, the Deterministic Model of COVID-19 is formulated. In the third Section, the model was analysed in terms of its positivity, invariant region, boundedness, existence, and uniqueness. In the fourth and fifth Sections, the Equilibrium Points and the Basic Reproduction Number were determined, and the stability of the Equilibrium Points was discussed, respectively. Also, the numerical results and discussion are mentioned in the sixth Section. Finally, the conclusions are given in Section 7.

This section presents positivity of the solution, the invariant region, and the existence and uniqueness of the solution for Model (1), by applying the technique used by Wameko et al. [25].

3.1 The solution's positivity

The starting condition of the model was assumed to be non-negative, and it is shown that the solution of the model is also positive.

Theorem 1. Solutions of the mathematical Model (1) together with the non-negative preliminary condition $P(0)=P_0, S(0)=S_0, I(0)=I_0$, and $T(0)=T_0$ are always positive; or the classes of model $P(t), S(t), I(t)$, and $T(t)$ are positive for all $t$.

Proof. Positivity of $P(t)$: Consider the first equation given by:

$\frac{d P}{d t}=\varphi \Lambda+\alpha T-(\delta+\beta) P$

After avoiding the positive parts $\varphi \Lambda+\alpha T$ from this equation, reduces to an inequality $\frac{d P}{d t} \geq$ $-(\delta+\beta) P$, and can be described in separable variables form as $\frac{d P}{P} \geq-(\delta+\beta) d t$. Following the integration and application of the initial condition $P(0)=P_0$, the solution can be acquired to be $P(t) \geq$ $P_0 e^{-(\delta+\beta) t}$. Now, since $P_0>0$ and $e^{-(\delta+\beta) t}>0$, it is simple to deduce that:

$P(t) \geq 0, \quad \forall t$.

Positivity of $S(t)$: Consider the second equation given by:

$\frac{d S}{d t}=(1-\varphi) \Lambda+\delta P+\eta T-(\varepsilon+\beta) S$

After avoiding the positive parts $(1-\varphi) \Lambda+\delta P+\eta T$ from this equation, reduces to an inequality $\frac{d S}{d t} \geq-(\varepsilon+\beta) S$ and can be described in separable variables form as $\frac{d S}{S} \geq-(\varepsilon+\beta) d t$. Following the integration and application of the initial condition $S(0)=S_0$, the solution can be acquired to be $S(t) \geq$ $S_0 e^{-(\varepsilon+\beta) t}$. Now, since $S_0>0$ and $e^{-(\varepsilon+\beta) t}>0$, it is simple to deduce that:

$S(t) \geq 0, \quad \forall t$.

Positivity of $I(t)$: By using the same method on the third equation, it can be obtained that:

$I(t) \geq 0, \quad \forall t$.

Positivity of $T(t)$: By using the same method on the fourth equation, it is possible to get that:

$T(t) \geq 0, \quad \forall t$.

We have shown a positive of the solution to the Model (1).

3.2 The invariant region

Theorem 2. The positive solutions of the system of model Eq. (1) are bounded; or, equivalently, the model variables $P(t), S(t), I(t)$, and $T(t)$ are bounded for all $t$. The Invariant Region in which the model solution is bounded has been obtained.

Proof. The entire population $\mathbb{N}$ in the model is given by:

$\mathbb{N}=P+S+I+T$

By deriving $\mathbb{N}$ in both directions with respect to $t$, which leads to

$\frac{d \mathbb{N}}{d t}=\frac{d P}{d t}+\frac{d S}{d t}+\frac{d I}{d t}+\frac{d T}{d t}$      (2)

By replacing model Eqs. (1) with (2), the following is obtained:

$\frac{d \mathbb{N}}{d t}=\Lambda-\beta \mathbb{N}-\gamma I$     (3)

If there are no COVID-19 related death i.e., $\gamma=0$, Eq. (3) becomes

$\frac{d \mathbb{N}}{d t}=\Lambda-\beta \mathbb{N} \Rightarrow \frac{1}{\Lambda-\beta \mathbb{N}} d \mathbb{N} \leq d t$,

and by integration the two sides, the following result is obtained:

$-\frac{1}{\beta} \ln (\Lambda-\beta \mathbb{N}) \leq t+c_1$,

where, $c_1$ is integration constant.

By using antilogarithms and inserting arbitrary constants $c_2=-\beta c_1$ and $B=e^{c_2}$. The generic solution may be written as follows:

$\Lambda-\beta \mathbb{N} \geq B e^{-\mu t}$     (4)

By putting preliminary conditions $\mathbb{N}(0)=\mathbb{N}_0$ in Eq. (4), the next step is obtained:

$B=\Lambda-\beta \mathbb{N}_0$      (5)

By substituting Eqs. (5) for (4), the result is the following:

$\Lambda-\beta \mathbb{N} \geq\left(\Lambda-\beta \mathbb{N}_0\right) e^{-\beta t} \Rightarrow \mathbb{N} \leq \frac{\Lambda}{\beta}-\frac{\left(\Lambda-\beta \mathbb{N}_0\right) e^{-\beta t}}{\beta}$      (6)

Accordingly, $t \rightarrow \infty$ in Eq. (6), the total population size is $\mathbb{N} \rightarrow \frac{\Lambda}{\beta}$, implying that $0 \leq \mathbb{N} \leq \frac{\Lambda}{\beta}$.

As a result, the model's viable solution set enters and remains in the region $\gamma=\{(P, S, I, T) \in$ $\left.\mathbb{R}_{+}^4: \mathbb{N} \leq \frac{\Lambda}{\beta}\right\}$

Therefore, the Model (1) is sound from an epidemiological and mathematical standpoint.

3.3 The solution's existence and uniqueness

Lemma 1. The solutions of the mathematical Model (1) together with the non-negative preliminary conditions $P(0)>0, S(0)>0, I(0)>0$, and $T(0)>0$ exist in $\mathbb{R}_{+}^4$, i.e., the model's solution for $P(t), S(t), I(t)$, and $T(t)$ exist for all $t$ and will continue to exist in $\mathbb{R}_{+}^4$.

Proof. The expressions on the right-hand side of Eq. (1) can be represented as follows:

$\begin{gathered}g_1(P, S, I, T)=\varphi \Lambda+\alpha T-(\delta+\beta) P \\ g_2(P, S, I, T)=(1-\varphi) \Lambda+\delta P+\eta T-(\varepsilon+\beta) S, \\ g_3(P, S, I, T)=\varepsilon S-(\gamma+\beta+\omega) I, \text { and } \\ g_4(P, S, I, T)=\omega I-(\eta+\alpha+\beta) T\end{gathered}$

Let $Y$ represent the region by Derrick and Grossman Theorem, $Y=\left\{(P, S, I, T) \in \mathbb{R}_{+}^4: \mathbb{N} \leq \frac{\Lambda}{\beta}\right\}$.

If $\frac{\partial g_i}{\partial y_j}, i, j=1,2,3,4$ are bounded and continuous in $Y$, then there is only one solution to the mathematical Model (1).

Also here, $y_1=P, y_2=S, y_3=I$ and $y_4=T$. The continuity and the boundedness are confirmed as follows:

In consequence, all of the partial derivatives $\frac{\partial g_i}{\partial y_j} ; i, j=1,2,3,4$ exist, bounded and continuous in $Y$.

Subsequently, the model's solution (1) is existing and unique, according to Derrick and Grossman Theorem.

In the coming section, the Equilibrium Points and Basic Reproduction Number $\mathcal{R}_0$ for Model (1) are demonstrated.

5.1 Stability of Disease-Free Equilibrium

5.1.1 Local stability of the Disease-Free Equilibrium

Proposition 1. If $\mathcal{R}_0<1$ Then the Disease-Free Equilibrium point is locally asymptotically- stable, and unstable if $\mathcal{R}_0>1$.

Proof. By following the approach outlined by Das [26], the proposition will be proved.

$J=\left[\begin{array}{cccc}-(\delta+\beta) & 0 & 0 & \alpha \\ \delta & -\left(\beta+\frac{(1-\vartheta) \pi \theta}{\mathbb{N}} I\right) & -\left(\frac{(1-\vartheta) \pi \theta}{\mathbb{N}}\right) S & \eta \\ 0 & \frac{(1-\vartheta) \pi \theta}{\mathbb{N}} I & \frac{(1-\vartheta) \pi \theta}{\mathbb{N}} S-(\gamma+\omega+\beta) & 0 \\ 0 & 0 & \omega & -(\eta+\alpha+\beta)\end{array}\right]$     (15)

$J_{\mathbb{E}_0}=\left[\begin{array}{cccc}-(\delta+\beta) & 0 & 0 & \alpha \\ \gamma & -\beta & -\frac{\pi \theta(1-\vartheta)(\delta+\beta-\varphi \beta)}{(\delta+\beta)} & \eta \\ 0 & 0 & \frac{\pi \theta(1-\vartheta)(\delta+\beta-\varphi \beta)}{(\delta+\beta)}-(\gamma+\beta+\omega) & 0 \\ 0 & 0& \omega & -(\eta+\alpha+\beta)\end{array}\right]$     (16)

The Jacobean Matrix is calculated at the DFE. $\mathbb{E}_0$, and the result is: $\operatorname{det}\left(J_{\mathbb{E}_0}-\lambda I\right)=0$, where $\lambda$ is eigenvalue, and $I=\left[\begin{array}{llll}1 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1\end{array}\right]$, the next step is obtained:

$\left|\begin{array}{cccc}-(\delta+\beta)-\lambda & 0 & 0 & \alpha \\ \delta & -\beta-\lambda & -\frac{\pi \theta(1-\vartheta)(\delta+\beta-\varphi \beta)}{(\delta+\beta)} & \eta \\ 0 & 0 & \frac{\pi \theta(1-\vartheta)(\delta+\beta-\varphi \beta)}{(\delta+\beta)}-(\gamma+\beta+\omega)-\lambda & 0 \\ 0 & 0 & \omega & -(\eta+\alpha+\beta)-\lambda\end{array}\right|=0$

$\Rightarrow(-(\beta+\delta)-\lambda)(-\beta-\lambda)\left(\frac{\pi \theta(1-\vartheta)(\delta+\beta-\varphi \beta)}{(\delta+\beta)}-(\beta+\omega+\gamma)-\lambda\right)(-(\beta+\eta+\alpha)-\lambda)=0$

Thus, the eigenvalues are $\lambda_1=-(\beta+\delta), \lambda_2=-\beta, \lambda_3=\frac{\pi \theta(1-\vartheta)(\delta+\beta-\varphi \beta)}{(\delta+\beta)}-(\beta+\omega+\gamma)$, and $\lambda_4=$ $-(\beta+\eta+\alpha)$

Hence, $\lambda_1, \lambda_2$, and $\lambda_4$ are negative, and $\lambda_3$ is negative if $\mathcal{R}_0<1$.

As a result, the Disease-Free Equilibrium point is locally asymptotically stable.

5.1.2 Global stability of the Disease-Free Equilibrium

This section utilises the methodology outlined by Kamgang and Sallet [27] to analyse the global stability of the DFE.

The equations given in Model (1) are expressed in the following manner:

$\begin{gathered}\frac{d X_s}{d t}=M\left(X_s-X_{D F E, s}\right)+M_1 X_i, \\ \frac{d X_i}{d t}=M_2 X_i .\end{gathered}$

where,

$M=\left[\begin{array}{ll}\frac{\partial g_1}{\partial P} & \frac{\partial g_1}{\partial S} \\ \frac{\partial g_2}{\partial P} & \frac{\partial g_2}{\partial S}\end{array}\right], M_1=\left[\begin{array}{ll}\frac{\partial g_1}{\partial I} & \frac{\partial g_1}{\partial T} \\ \frac{\partial g_2}{\partial I} & \frac{\partial g_2}{\partial T}\end{array}\right]$, and $M_2=\left[\begin{array}{ll}\frac{\partial g_3}{\partial I} & \frac{\partial g_3}{\partial T} \\ \frac{\partial g_4}{\partial I} & \frac{\partial g_4}{\partial T}\end{array}\right]$,

such that

$\begin{gathered}g_1(P, S, I, T)=\varphi \Lambda+\alpha T-(\delta+\beta) P, \\ g_2(P, S, I, T)=(1-\varphi) \Lambda+\delta P+\eta T-(\varepsilon+\beta) S, \\ g_3(P, S, I, T)=\varepsilon S-(\gamma+\beta+\omega) I, \\ g_4(P, S, I, T)=\omega I-(\eta+\alpha+\beta) T .\end{gathered}$

The vector corresponding to the transmitting classes is $X_i$, and the vector corresponding the nontransmitting class is $X_s$.

If $M$ has negative. eigenvalues, and $M_2$ is a Metzler Matrix, meaning that its off-diagonal elements do not have negative values, then the DFE is globally asymptotically stable.

For the Model Eq. (1), we have $X_s=(P, S)^t$, and $X_i=(I, T)^t$, where a transposition of the matrix is indicated by the superscript $t$.

It is necessary to confirm that $M_2$ is a Metzler Matrix, and that the matrix $M$ for non-transmitting classes indeed has negative eigenvalues.

The next is obtained from the Model's Eq. (1) for non-transmitting classes:

$M=\left[\begin{array}{cc}-(\delta+\beta) & 0 \\ \delta & -(\varepsilon+\beta)\end{array}\right]$

The eigenvalues are obtained from the matrix $M$ above, $\lambda_1=-(\delta+\beta)$, and $\lambda_2=-(\varepsilon+\beta)$, all of them are real and negative.

Consequently, the system must be $\frac{d X_s}{d t}=M\left(X_s-X_{D F E, s}\right)+M_1 X_i$ at DFE is locally and globally asymptotically stable if $M_2$ is Metzler Matrix.

Rearranging in the right way yields:

$\begin{gathered}M_2=\left[\begin{array}{cc}\left(\frac{\pi \theta(1-\vartheta)}{\mathbb{N}} S-(\gamma+\beta+\omega)\right. & 0 \\ \omega & -(\eta+\alpha+\beta)\end{array}\right], \\ M_1=\left[\begin{array}{cc}0 & \alpha \\ -\left(\frac{\pi \theta(1-\vartheta)}{\mathbb{N}} S\right. & \eta\end{array}\right], \\ X_s-X_{D F E, s}=\left(\begin{array}{c}P-\frac{\varphi \Lambda}{(\delta+\beta)} \\ S-\frac{(\delta+\beta-\varphi \beta) \Lambda}{\beta(\delta+\beta)}\end{array}\right)\end{gathered}$

$M_2$ is a Metzler Matrix because it has non-negative off-diagonal elements.

Consequently, the DFE is asymptotically globally stable.

5.2 Stability of the Endemic Equilibrium

5.2.1 Global stability of the Endemic Equilibrium

Theorem 3. The Endemic. Equilibrium points of the Model (1) are globally asymptotically stable in $Y$ if $\mathcal{R}_0>1$.

Proof. This theorem can be proved using the method described by Tilahun et al. [28].

The Lyapunov function is utilized as follows:

$\mathcal{L}\left(P^*, S^*, I^*, T^*\right)=\left(P-P^*-P^* \ln (P)\right)+\left(S-S^*-S^* \ln (S)\right)+\left(I-I^*-I^* \ln (I)\right)+\left(T-T^*-T^* \ln (T)\right)$

By derivative $\mathcal{L}$ with respect to $t$, the next step is obtained:

$\frac{d \mathcal{L}}{d t}=\left(1-\frac{P^*}{P}\right) \frac{d P}{d t}+\left(1-\frac{S^*}{S}\right) \frac{d S}{d t}+\left(1-\frac{I^*}{I}\right) \frac{d I}{d t}+\left(1-\frac{T^*}{T}\right) \frac{d T}{d t}$

By substituting Model (1) into the above step, the following is obtained:

$\begin{aligned} & \frac{d \mathcal{L}}{d t}=\left(1-\frac{P^*}{P}\right)(\varphi \Lambda+\alpha T-(\delta+\beta) P)+\left(1-\frac{S^*}{S}\right)((1-\varphi) \Lambda+\delta P+\eta T-(\varepsilon+\beta) S) \\ & +\left(1-\frac{I^*}{I}\right)(\varepsilon S-(\gamma+\beta+\omega) I)+\left(1-\frac{T^*}{T}\right)(\omega I-(\eta+\alpha+\beta) T) \\ & \frac{d \mathcal{L}}{d t}=\varphi \Lambda+\alpha T-(\delta+\beta) P-\frac{P^*}{P}(\varphi \Lambda+\alpha T)+(\delta+\beta) P^*+(1-\varphi) \Lambda+\delta P+\eta T-(\varepsilon+\beta) S \\ & -\frac{S^*}{S}((1-\varphi) \Lambda+\delta P+\eta T)+(\varepsilon+\beta) S^*+\varepsilon S-(\gamma+\beta+\omega) I-\frac{I^*}{I} \varepsilon S+(\gamma+\beta+\omega) I^*+\omega I \\ & -(\eta+\alpha+\beta) T-\frac{T^*}{T} \omega I+(\eta+\alpha+\beta) T^* \\ & =\Lambda+(\delta+\beta) P^*+(\varepsilon+\beta) S^*+(\gamma+\beta+\omega) I^*+(\eta+\alpha+\beta) T^* \\ & -\left[\frac{P^*}{P}(\varphi \Lambda+\alpha T)+\frac{S^*}{S}((1-\varphi) \Lambda+\delta P+\eta T)+\frac{I^*}{I} \varepsilon S+\frac{T^*}{T} \omega I+\beta P+\beta S+(\gamma+\beta) I+\beta T\right] \\ & \end{aligned}$     (17)

Let

$\begin{gathered}\mathcal{L}_1=\Lambda+(\delta+\beta) P^*+(\varepsilon+\beta) S^*+(\gamma+\beta+\omega) I^*+(\eta+\alpha+\beta) T^* \\ \mathcal{L}_2=\frac{P^*}{P}(\varphi \Lambda+\alpha T)+\frac{S^*}{S}((1-\varphi) \Lambda+\delta P+\eta T)+\frac{I^*}{I} \varepsilon S+\frac{T^*}{T} \omega I+\beta P+\beta S+(\gamma+\beta) I+\beta T\end{gathered}$

Then Eq. (17) is written as follows:

$\frac{d \mathcal{L}}{d t}=\mathcal{L}_1-\mathcal{L}_2$     (18)

From Eq. (18) $\frac{d \mathcal{L}}{d t} \leq 0$, when $\mathcal{L}_1 \leq \mathcal{L}_2$, and $\frac{d \mathcal{L}}{d t}=0$ if and only if $P=P^*, S=S^*, I=I^*$, and $T=T^*$.

The largest compact invariant set in $Y=\left\{(P, S, I, T) \in \mathbb{R}_{+}^4: \frac{d \mathcal{L}}{d t}=0\right\}$ is the singleton of EE.

Therefore, the Endemic Equilibrium point is globally asymptotically stable in $Y$.

In the subsequent section, we will introduce numerical results and discussion for Model (1).

In recent times, the COVID-19 pandemic has arisen as a substantial concern within the realm of the general population's health. The PSITPS Model for the COVID-19 epidemic has been. studied. The existence. and uniqueness. of the. solution, as well as its positivity and boundedness have been proven for the proposal Model. The determination of the Basic. Reproduction. Number $\mathcal{R}_0$ is crucial in assessing the potential for disease eradication or spread within a population. The study conducted a stability. analysis, which revealed that the Disease-Free Equilibrium exhibits local asymptotic stability under the condition that $\mathcal{R}_0<1$. Conversely, the Endemic. Equilibrium point is established when $\mathcal{R}_0$ exceeds 1 .

The parameter values for the model were determined through the utilization of actual data obtained from the Public Health Department and the Statistics Division of the Training and Human Development Center of the Najaf Health Directorate. The model studied in this work was run through numerical simulations using the RK4 method in the MATLAB program. The findings of this study suggest that the implementation of protective measures, such as vaccination, significantly contributes to the mitigation of COVID-19 transmission in the region of Najaf for a prolonged duration. Furthermore, it is imperative to implement measures that decrease the contact rate while simultaneously increasing the vaccination rate to effectively manage the disease.