The Analytic Solutions of the Fractional-Order Model for the Spatial Epidemiology of the COVID-19 Infection

Abstract

This paper provides a mathematical fractional-order model that accounts for the mindset of patients in the transmission of COVID-19 disease, the continuous inflow of foreigners into the country, immunization of population subjects, and temporary loss of immunity by recovered individuals. The analytic solutions, which are given as series solutions, are derived using the fractional power series method (FPSM) and the residual power series method (RPSM). In comparison, the series solution for the number of susceptible members, using the FPSM, is proportional to the series solution, using the RPSM for the first two terms, with a proportional constant of $\psi \Gamma n\alpha +1$, where $\psi$ is the natural birth rate of the baby into the susceptible population, $\Gamma$ is the gamma function, $n$ is the $n$th term of the series, and $\alpha$ is the fractional order as the initial number of susceptible individuals approaches the population size of Ghana. However, the variation in the two series solutions of the number of members who are susceptible to the COVID-19 disease begins at the third term and continues through the remaining terms. This is brought on by the nonlinear function present in the equation for the susceptible subgroup. The similar finding is made in the series solution of the number of exposed individuals. The series solutions for the number of deviant people, the number of nondeviant people, the number of people quarantined, and the number of people recovered using the FPSM are unquestionably almost identical to the series solutions for same subgroups using the RPSM, with the exception that these series solutions have initial conditions of the subgroup of the population size. It is observed that, in this paper, the series solutions of the nonlinear system of fractional partial differential equations (PDEs) provided by the RPSM are more in line with the field data than the series solutions provided by the FPSM.

This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1. Introduction

The development of a mathematical model for understanding and unravelling the underlying mechanisms of the epidemiology of the COVID-19 disease has garnered interest from public health systems to academia in several different countries; the majority of these models focus on epidemics of the disease progression from one person to another person, as described by [1–3]. However, the COVID-19 disease originated in Wuhan, China, and geographically spread to other parts of the world as a pandemic disease [4]. In this sense, COVID-19 epidemiology is more appropriately classified as a pandemic disease than an epidemic disease. The primary factor in the global transmission of the COVID-19 disease is the movement of exposed or infected individuals, who may or may not have the aim of coming into contact with the vulnerable individuals in the host country. The spatial spread of the disease in the various countries was accounted for in the mathematical models developed by [5, 6] by taking into consideration the diffusing susceptible individuals, exposed individuals, and infected individuals. Despite this, these models do not account for the vaccinations that the individuals of the population received.

There are currently treatments available that are given to people all around the world regardless of their health conditions, such as the Johnson and Johnson vaccine and the AstraZeneca vaccine. Although the usefulness of these vaccines has been scientifically demonstrated, these immunizations do lose some of their efficacy over time. There is no assurance that a person receiving the COVID-19 vaccination will be protected from getting the disease upon contact with a person with the SARS virus. For example, see authors in [7]. This observation makes it unclear which individuals are completely unprotected from the disease and which persons are temporarily protected for a short period of time following immunization. Only a few researchers have used vaccinated subjects in their models without the inclusion of the spatial transmission of the disease. For example, have a look at the authors in [8–10]. All of these epidemiological models account for people moving from one subgroup of the population to another subgroup of the population. Although they captured vaccinated persons, they did not incorporate the diffusing individuals who brought the COVID-19 disease into their respective countries. Also, for the mathematical models on the control of transmission of COVID-19 disease, see [11]. The findings of these researchers, however, are not all inclusive since they neglected to take into account an important observation of the progression of patient through the disease. Evidence from numerous countries has demonstrated that infected individuals (patients) either plan or do not intend to transmit the COVID-19 virus to susceptible individuals [4]. In developing a mathematical model to describe the epidemiology of the COVID-19 infection, the mindset of the spreaders was not captured in their models. Additionally, statistics from different countries have revealed that whether a person takes medicine to treat COVID-19 or not, they still run the risk of getting the illness again if they come into contact with an infected person. Thus, the recovery from the disease is for a short period of time (see [12]). When creating a mathematical model to describe the epidemiology of COVID-19, all these issues were not taken into account.

The type of mathematical tools a researcher uses to arrive at his or her conclusion(s) ultimately determines the success of any mathematical analysis. Since the beginning of the COVID-19 outbreak in China till now, researchers have mainly relied heavily on either the use of the qualitative method, the quantitative method, or both. These methods have significantly more drawbacks than advantages. A numerical scheme is an example of a quantitative method that always approximates the exact solution of the differential equation with some level of precision. This quantitative method yields intolerable inaccuracies; in the worst situation, its solution diverges from the exact solution of a differential equation. As usual, even if the solution suggested by the numerical scheme exists, one must perform a number of iterations before reaching the desired solution. The qualitative method narrows down the information contained in the solution of the differential equation. The domain elements of the function that describes the epidemiology of COVID-19 infection are revealed by this method of investigation on a microscopic level. In light of this, this method can only produce fixed-point solutions of differential equations. The vast majority of nonstationary points are uncovered by this method. More crucially, neither a quantitative nor a qualitative method provides the function that describes the theoretical foundation for describing the epidemiology of COVID-19 disease.

Recently, it has been discovered that the integer differential equations suffer from several shortcomings when compared to the differential equations of fractional order. The fractional differential equation has memory and heredity properties because of its nonlocal property for describing COVID-19 pandemics. Any solution to the system of PDEs, regardless of order, may be easily obtained using fractional ordering. The theory of controls, infectious diseases, growth of tumours, and feedback systems are examples of applied scientific problems where the differential equations of fractional order have proven to be effective models. For example, see authors in [13] who applied the fractal fraction Adams-Bashforth method to search for the solution of fractal-fractional susceptible-infective-recovered model. Another numerical approach for solving systems of differential equations that are both linear and nonlinear is the Pade approximation method. High-order approximations are necessary when using this method. More crucially, given a nonlinear system of PDEs, there is no systematic procedure in selecting the parameters in the Pade approximation method[14]. Since Mittag-Leffler functions or their derivatives make up the majority of the solutions to the system of fractional differential equations, rigorous mathematics is necessary to solve these equations. One of the methods for solving system of fractional differential equations is the RPSM which was first observed by [15] for solving fuzzy differential equation. With this approach, a power series is assumed to exist for the system of ODEs, and the coefficients of the power series are used to create a recurrence equation. When the residual coefficients, from the power series, are equal to zero, an algebraic system of equations results, from which the values of the series solution of the unknown coefficients can be deduced. Nevertheless, while solving fractional-order PDE, say in two variables, this method assumes that one of the independent variables has a representation in a fractional power series, and the second independent variable is handled as a coefficient variable, which is roughly derived from the variation in the given fractional-order PDE based on the initial or boundary condition. For example, see authors in [16–19]. The same method was used by [20] to solve nonlinear fractional-order PDEs. In [21], the authors used the Atangana-Baleanu fractional derivative to obtain asymptotic interval approximation solutions to the fractional differential equation under various conditions. A nonlinear system of stiff fractional-order PDEs and the nonlinear system of fractional PDEs have not been solved using the RPSM. The kind of nonlinearity in a fractional PDE largely depends on the functional space which contains the solution of fractional differential PDE. The FPSM is another intriguing method which was first observed by [22]. The authors in [23, 24] applied this method to obtain the solutions of fractional PDEs. It is challenging to find analytic solutions to a nonlinear system of fractional-order partial differential equations. Additionally, researchers from all over the world have not observed a comparison of the series solutions utilizing both the RPSM and the FPSM. More importantly, no information regarding comparing the series solutions obtained by these methods with field data is provided in the literature. The series solution (analytic) method of the nonlinear system of fractional-order partial differential equations has a solution, is the most dependable and efficient method as compared to both the qualitative and the quantitative methods.

In this paper, the infected group of the SEIQR model is further divided into two subgroups: the deviant infected subgroup and nondeviant subgroup of the population in the classical susceptible-exposed-infected-quarantined-recovered model with diffusion terms. Thus, the susceptible-exposed-deviant infected-nondeviant infected-quarantined-recovered (SEII${}^{\ast }$QR) model with diffusion terms, and vaccinated susceptible term is developed. In addition, the fractional form of this model is provided herein. Moreover, both the FPSM and the RPSM are used to obtain the series solution (analytic) of the nonlinear system of fractional PDEs. The solutions that are yielded by these two methods are compared with field data accounting for the robustness of the methods.

2. Fundamental Concept in Fractional Calculus

Definition 1.

A real function $ux,t,x\in I,t>0$ is said to be in the space $C_{\alpha }I\times {ℝ}^{+},\mu \in ℝ$, if there exist a real number $p>\alpha$ such that $ux,t=t^{p}fx,t$, where $fx,t\in CI\times {ℝ}^{+}$, and it is said to be in the space $C_{\alpha }^n$, if ${\partial }^n/\partial t^n\in C_{\alpha }$, $n\in \mathbb{N}$ (see [25]).

Definition 2.

For $n-1<\beta <n$, $n\in \mathbb{N}$. The Caputo fractional derivative operator of the order $\beta$ is define by (see [26]) $$\begin{array}{}\text{(1)}& D_a^{\alpha }ut=\frac{1}{\Gamma n-\alpha }{\int }_a^t{u}^n\xi {t-\xi }^{n-\alpha -1}d\xi ,t>0.\end{array}$$

Theorem 3.

The fractional power series (FPS) ${\sum }_{m=0}^{\infty }{a}_m{t-t_0}^{m\alpha }$:

(i) converges only for $t=t_0$, that is, the radius of convergence equal to zero

(ii) converges for all $t\ge t_0$, that is; the radius of convergence equal to $\infty$

(iii) converges for $t\in t_0,t_0+R$, for some positive real numbers $R$, and diverges for $t>t_0+R$. Here, $R$ is the radius of convergence for the FPS [27].

3. Main Results

In this section, for modelling the COVID-19 epidemiology in Ghana, a mathematical model that takes into account the mindset of the patients in spreading the COVID-19 disease, temporary loss of immunity by recoveries, and the continuous influx of foreigners entering the country with or without the disease is needed. Therein, the FPSM- and RPSM-based analytical solutions of the nonlinear system of fractional PDEs are presented as series solutions.

3.1. Model Description

Despite the fact that the COVID-19 vaccination is given to country residents by the Ministry of Health (MOH), neither the Johnson and Johnson nor the AstraZeneca vaccine is anticipated to provide COVID-19 patients with a lifetime of immunity against the illness.

In Figure 1, the population size of Ghana is split into six distinct subgroups namely: the susceptible subgroup, $Sx,t$; exposed subgroup, $Ex,t$; deviant infected subgroup, $Ix,t$; nondeviant infected subgroup, $I^{\ast }x,t$; quarantined subgroup, $Qt$; and the recovered subgroup, $Rt$. A susceptible person is any member of the population who is capable of catching the SARS virus from an infected COVID-19 patient. An exposed person is someone who has caught the SARS virus, but for a brief while, he or she is unable to pass it on to a susceptible person. The waiting period is therefore in effect for this person. The deviant infected person is a patient who has chosen to purposefully spread the SARS virus to susceptible family members on the grounds that since they already have the illness, they must also experience the COVID-19 disease-related consequences. The nondeviant individual, on the other hand, is an infected person who does not willingly spread the SARS virus to susceptible family members or friends because they do not want them to get the COVID-19 illness. The person under quarantine is a COVID-19 patient who was deviant infected person or nondeviant infected person, or an exposed person, whose movements are restricted in a specific place for an extended length of time. The recovered person is the person who either recovers naturally or receives treatment at the hospital for a period of time. After being exposed to the COVID-19 infection, this person is still at a significant risk of reacquiring the SARS virus. As a consequence of this, the immune system is unable to recuperate and is therefore vulnerable to losing its immunity. To take into account the continuous inflow of foreigners into the country, the subgroups of susceptible, exposed, deviant infected person, and nondeviant infected person depend on the distance, $x$, as well as the passage of time.

[figure(s) omitted; refer to PDF]

Additionally, $\psi$ stands for the natural birth rate, $\beta$ for transmission rate, and $\mu$ for natural death rate, $\nu$ is the rate for quarantining exposed individuals, and ${\alpha }_1$ is the rate at which an exposed person intends to transmit the SARS virus to susceptible members. The rate at which an exposed person has no intention to infect a susceptible member with the SARS virus is ${\alpha }_2$. The rates at which the exposed person, the deviant infected person, and the nondeviant infected person are quarantined, respectively, are $\nu$, ${\omega }_1$, and ${\omega }_2$. The disease-induced death rates from the subgroups of deviants, nondeviants, and confined individuals are ${\delta }_1$, ${\delta }_2$, and ${\delta }_3$, respectively. The $\varphi$, ${\gamma }_1$, and ${\gamma }_2$ are the rates of recoveries from COVID-19 disease by the quarantined, the deviant, and the nondeviant individuals, respectively. The $\eta$ is the rate at which susceptible moves to the recovered compartment after receiving a vaccine, and $\kappa$ is the rate at which recovered person loses their immunity and becomes susceptible again. The natural death rate from each subgroup of the population is denoted by $\mu$. Due to the fact that the model takes into account diffusion of the foreigners into the susceptible, exposed, deviant, and nondeviant subgroups, the rates of diffusion into the corresponding compartments are specified as follows: $D_1$ represents the rate of diffusion into the susceptible compartment, $D_2$ represents the rate of diffusion into the exposed compartment, $D_3$ represents the rate of diffusion into the deviant subgroup, and $D_4$ represents the rate of diffusion into the nondeviant subgroup.

Based on above facts, the following nonlinear system of fractional PDEs is obtained for describing the epidemiology of COVID-19 in Ghana. $$\begin{array}{}\text{(2)}& \frac{{\partial }^{\alpha }S}{\partial t^{\alpha }}=\psi N+\kappa R-\frac{\beta SI+I^{\ast }+Q}{N}-\mu +\eta S+D_1\frac{{\partial }^{2\alpha }S}{\partial x^{2\alpha }},\text{(3)}& \frac{{\partial }^{\alpha }E}{\partial t^{\alpha }}=\frac{\beta SI+I^{\ast }+Q}{N}-{\alpha }_1+{\alpha }_2+\mu +vE+D_2\frac{{\partial }^{2\alpha }E}{\partial x^{2\alpha }},\text{(4)}& \frac{{\partial }^{\alpha }I}{\partial t^{\alpha }}={\alpha }_{1}E-{\gamma }_1+\mu +{\delta }_1+{\omega }_{1}I+D_3\frac{{\partial }^{2\alpha }I}{\partial x^{2\alpha }},\text{(5)}& \frac{{\partial }^{\alpha }{I}^{\ast }}{\partial t^{\alpha }}={\alpha }_{2}E-{\gamma }_2+\mu +{\delta }_2+{\omega }_2{I}^{\ast }+D_4\frac{{\partial }^{2\alpha }{I}^{\ast }}{\partial x^{2\alpha }},\text{(6)}& \frac{d^{\alpha }Q}{d{t}^{\alpha }}=vE+{\omega }_{1}I+{\omega }_2{I}^{\ast }-\mu +\varphi +{\delta }_{3}E,\text{(7)}& \frac{{\partial }^{\alpha }R}{\partial t^{\alpha }}=\varphi Q+{\gamma }_{1}I+{\gamma }_2{I}^{\ast }+\eta S-\mu +\kappa R,\end{array}$$together with the initial conditions $$\begin{array}{c}\text{(8)}& Sx,0=S_{0}x,Ex,0=E_{0}x,\\ Ix,0=I_{0}x,\\ I^{\ast }x,0=I_0^{\ast }x,\\ Q0=0,\\ R0=0.\end{array}$$

3.2. Analytic Solutions of the System of Fractional Partial Differential Equations Using the Fractional Power Series Method

In this section, series solutions of the system of equations (2)–(7) together with the initial conditions in equation (8) are obtained in Hilbert space using the FPSM. In obtaining each solution of the system of equations (2)–(7) together with initial conditions, it is assumed that the unknown function defining the equation is in series form which converges to a known function. In addition, the proof of the existence of these series solution as well as its uniqueness is provided here.

Setting $$\begin{array}{}\text{(9)}& Sx,t=\sum _{k=0}^{\infty }{S}_{k}x{t}^{k\alpha },\text{(10)}& Ex,t=\sum _{k=0}^{\infty }{E}_{k}x{t}^{k\alpha },\text{(11)}& Ix,t=\sum _{k=0}^{\infty }{I}_{k}x{t}^{k\alpha },\text{(12)}& I^{\ast }x,t=\sum _{k=0}^{\infty }{I}_k^{\ast }x{t}^{k\alpha },\text{(13)}& Qx,t=\sum _{k=0}^{\infty }{Q}_{k}x{t}^{k\alpha },\text{(14)}& Rx,t=\sum _{k=0}^{\infty }{R}_{k}x{t}^{k\alpha },\end{array}$$we can see that $$\begin{array}{c}\text{(15)}& D_t^{\alpha }Sx,t=\sum _{k=1}^{\infty }\frac{\Gamma k\alpha +1}{\Gamma k-1\alpha +1}{S}_{k}x{t}^{k-1\alpha },\frac{{\partial }^{2\alpha }S}{\partial x^{2\alpha }}=\sum _{k=0}^{\infty }\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{S}_{k}x{t}^{k\alpha }.\end{array}$$

Substituting equations (8), (9), (10), (11), (12), (13), and (14) into equation (2) yields $$\begin{array}{}\text{(16)}& \sum _{k=1}^{\infty }{S}_{k}x\frac{\Gamma k\alpha +1}{\Gamma k-1\alpha +1}{t}^{k-1\alpha }=\psi N+\kappa \sum _{k=0}^{\infty }{R}_{k}x{t}^{k\alpha }-\frac{\beta }{N}\sum _{k=0}^{\infty }{S}_{k}x{t}^{k\alpha }\sum _{k=0}^{\infty }{I}_{k}x{t}^{k\alpha }+\sum _{k=0}^{\infty }{I}_k^{\ast }x{t}^{k\alpha }+\sum _{k=0}^{\infty }{Q}_{k}x{t}^{k\alpha }-\mu +\eta \sum _{k=0}^{\infty }{S}_{k}x{t}^{k\alpha }+D_1\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}\sum _{k=0}^{\infty }{S}_{k}x{t}^{k\alpha }.\end{array}$$

Comparing the powers of $t^0$, we have $$\begin{array}{}\text{(17)}& S_{1}x=\frac{1}{\Gamma \alpha +1}\psi N+\kappa R_{0}x-\frac{\beta }{N}{S}_{0}x{I}_{0}x+I_0^{\ast }x+Q_{0}x-\mu +\eta S_{0}x+D_1\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{S}_{0}x.\end{array}$$

Similarly, we observe the following results. For $S_{2}x,S_{3}x,\cdots ,S_{n}x$, we have $$\begin{array}{c}\text{(18)}& S_{2}x=\frac{\Gamma \alpha +1}{\Gamma 2\alpha +1}\psi N+\kappa R_{1}x-\frac{\beta }{N}{S}_{0}x{I}_{1}x+S_{0}x{I}_1^{\ast }x+S_{0}x{Q}_{1}x+S_{1}x{I}_{0}x+S_{1}x{I}_0^{\ast }x+S_{1}x{Q}_{0}x-\mu +\eta S_1+D_1\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{S}_1,S_{3}x=\frac{\Gamma 2\alpha +1}{\Gamma 3\alpha +1}\psi N+\kappa R_{2}x-\frac{\beta }{N}{S}_{0}x{I}_{2}x+S_{0}x{I}_2^{\ast }x+S_{0}x{Q}_{2}x+S_{1}x{I}_{1}x+S_{1}x{I}_1^{\ast }x+S_{1}x{Q}_{1}x+S_{2}x{I}_{0}x+S_{2}x{I}_0^{\ast }x+S_{2}x{Q}_{0}x-\mu +\eta S_2+D_1\frac{{\partial }^{2\alpha }{S}_2}{\partial x^{2\alpha }},\\ ⋮\\ S_{n}x=\frac{\Gamma n-1\alpha +1}{\Gamma n\alpha +1}\psi N+\kappa R_{n-1}x-\frac{\beta }{N}\sum _{i=0}^{n-1}{S}_{i}x{I}_{n-1-i}x+S_{i}x{I}_{n-1-i}^{\ast }x+S_{i}x{Q}_{n-1-i}x-\mu +\eta S_{n-1}x+D_1\frac{{\partial }^{2\alpha }{S}_{n-1}x}{\partial x^{2\alpha }}.\end{array}$$

Following the similar procedure above, the following results are obtained for $E_{1}x,E_{2}x,E_{3}x,\cdots ,E_{n}x$: $$\begin{array}{c}\text{(19)}& E_{1}x=\frac{1}{\Gamma \alpha +1}\frac{\beta }{N}{S}_{0}x{I}_{0}x+I_0^{\ast }x+Q_{0}x-{\alpha }_1+{\alpha }_2+\mu +\nu E_0+D_2\frac{{\partial }^{2\alpha }{E}_{0}x}{\partial x^{2\alpha }},E_{2}x=\frac{\Gamma \alpha +1}{\Gamma 2\alpha +1}\frac{\beta }{N}{S}_{0}x{I}_{1}x+S_{0}x{I}_1^{\ast }x+S_{0}x{Q}_{1}x+I_{0}x{S}_{1}x+S_{1}x{I}_0^{\ast }x+S_{1}x{Q}_{0}x-{\alpha }_1+{\alpha }_2+\mu +\nu E_1+D_2\frac{{\partial }^{2\alpha }{E}_{1}x}{\partial x^{2\alpha }},\\ E_{3}x=\frac{\Gamma 2\alpha +1}{\Gamma 3\alpha +1}\frac{\beta }{N}{S}_{0}x{I}_{2}x+S_{0}x{I}_2^{\ast }x+S_{0}x{Q}_{2}x+S_{1}x{I}_{1}x+S_{1}x{I}_1^{\ast }x+S_{1}x{Q}_{1}x+S_{2}x{I}_{0}x+S_{2}x{I}_0^{\ast }x+S_{2}x{Q}_{0}x-{\alpha }_1+{\alpha }_2+\mu +\nu E_2+D_2\frac{{\partial }^{2\alpha }{E}_{2}x}{\partial x^{2\alpha }},\\ ⋮\\ E_{n}x=\frac{\Gamma n-1\alpha +1}{\Gamma n\alpha +1}\frac{\beta }{N}\sum _{k=0}^{n-1}{S}_{k}x{I}_{n-1-k}x+S_{k}x{I}_{n-1-k}^{\ast }x+S_{k}x{I}_{n-1-k}x-{\alpha }_1+{\alpha }_2+\mu +\nu E_{n-1}x+D_2\frac{{\partial }^{2\alpha }{E}_{n-1}xx}{\partial x^{2\alpha }}.\end{array}$$

Similarly, the series solutions for the number of the deviant infected people, the number of nondeviant infected people, the quarantined, and the number of recoveries are as follows: $$\begin{array}{c}\text{(20)}& I_{1}x=\frac{1}{\Gamma \alpha +1}{\alpha }_1{E}_{0}x-{\gamma }_1+\mu +{\delta }_1+{\omega }_1{I}_{0}x+D_3\frac{{\partial }^{2\alpha }{I}_{0}x}{\partial x^{2\alpha }},I_{2}x=\frac{\Gamma \alpha +1}{\Gamma 2\alpha +1}{\alpha }_1{E}_{1}x-{\gamma }_1+\mu +{\delta }_1+{\omega }_1{I}_{1}x+D_3\frac{{\partial }^{2\alpha }{I}_{1}x}{\partial x^{2\alpha }},\\ I_{3}x=\frac{\Gamma 2\alpha +1}{\Gamma 3\alpha +1}{\alpha }_1{E}_{2}x-{\gamma }_1+\mu +{\delta }_1+{\omega }_1{I}_{2}x+D_3\frac{{\partial }^{2\alpha }{I}_{2}x}{\partial x^{2\alpha }},\\ ⋮\\ I_{n}x=\frac{\Gamma n-1\alpha +1}{\Gamma n\alpha +1}{\alpha }_1{E}_{n-1}x-{\gamma }_1+\mu +{\delta }_1+{\omega }_1{I}_{n-1}x+D_3\frac{{\partial }^{2\alpha }{I}_{n-1}x}{\partial x^{2\alpha }},\\ \text{(21)}& I_1^{\ast }x=\frac{1}{\Gamma \alpha +1}{\alpha }_2{E}_{0}x-{\gamma }_2+\mu +{\delta }_2+{\omega }_2{I}_0^{\ast }x+D_4\frac{{\partial }^{2\alpha }{I}_0^{\ast }x}{\partial x^{2\alpha }},I_2^{\ast }x=\frac{\Gamma \alpha +1}{\Gamma 2\alpha +1}{\alpha }_2{E}_{1}x-{\gamma }_2+\mu +{\delta }_2+{\omega }_2{I}_1^{\ast }x+D_4\frac{{\partial }^{2\alpha }{I}_1^{\ast }x}{\partial x^{2\alpha }},\\ I_3^{\ast }x=\frac{\Gamma 2\alpha +1}{\Gamma 3\alpha +1}{\alpha }_2{E}_{2}x-{\gamma }_2+\mu +{\delta }_2+{\omega }_2{I}_2^{\ast }x+D_4\frac{{\partial }^{2\alpha }{I}_2^{\ast }x}{\partial x^{2\alpha }},\\ ⋮\\ I_n^{\ast }x=\frac{\Gamma n-1\alpha +1}{\Gamma n\alpha +1}{\alpha }_2{E}_{n-1}x-{\gamma }_2+\mu +{\delta }_2+{\omega }_2{I}_{n-1}^{\ast }x+D_4\frac{{\partial }^{2\alpha }{I}_{n-1}^{\ast }x}{\partial x^{2\alpha }},\\ \text{(22)}& Q_{1}x=\frac{1}{\Gamma \alpha +1}v{E}_{0}x+{\omega }_1{I}_{0}x+{\omega }_2{I}_0^{\ast }x-\mu +\varphi +{\delta }_3{Q}_{0}x,Q_{2}x=\frac{\Gamma \alpha +1}{\Gamma 2\alpha +1}v{E}_{1}x+{\omega }_1{I}_{1}x+{\omega }_2{I}_1^{\ast }x-\mu +\varphi +{\delta }_3{Q}_{1}x,\\ Q_{3}x=\frac{\Gamma 2\alpha +1}{\Gamma 3\alpha +1}v{E}_{2}x+{\omega }_1{I}_{2}x+{\omega }_2{I}_2^{\ast }x-\mu +\varphi +{\delta }_3{Q}_{2}x,\\ ⋮\\ Q_{n}x=\frac{\Gamma n-1\alpha +1}{\Gamma n\alpha +1}v{E}_{n-1}+{\omega }_1{I}_{n-1}x+{\omega }_2{I}_{n-1}^{\ast }x-\mu +\psi +{\delta }_3{Q}_{n-1}x,\\ \text{(23)}& R_{1}x=\frac{1}{\Gamma \alpha +1}\varphi Q_0+{\gamma }_1{I}_0+{\gamma }_2{I}_0^{\ast }+\eta S_0-\mu +\kappa R_0,R_{2}x=\frac{\Gamma \alpha +1}{\Gamma 2\alpha +1}\varphi Q_1+{\gamma }_1{I}_1+{\gamma }_2{I}_1^{\ast }+\eta S_1-\mu +\kappa R_1,\\ R_{3}x=\frac{\Gamma 2\alpha +1}{\Gamma 3\alpha +1}\varphi Q_2+{\gamma }_1{I}_2+{\gamma }_2{I}_2^{\ast }+\eta S_2-\mu +\kappa R_2,\\ ⋮\\ R_{n}x=\frac{\Gamma n-1\alpha +1}{\Gamma n\alpha +1}\varphi Q_{n-1}x+{\gamma }_1{I}_{n-1}x+{\gamma }_2{I}_{n-1}^{\ast }x+\eta S_{n-1}x-\mu +\kappa R_{n-1}x.\end{array}$$

The series solutions of the nonlinear system of fractional PDEs order of $n$th term are given by $$\begin{array}{}\text{(24)}& S_{n}x,t=\sum _{n=1}^{\infty }\frac{\Gamma n-1\alpha +1}{\Gamma n\alpha +1}\psi N+\sum _{n=1}^{\infty }\kappa R_{n-1}x-\frac{\beta }{N}\sum _{k=0}^{n-1}{S}_{k}x{I}_{n-1-k}x+S_{k}x{I}_{n-1-k}^{\ast }x+S_{k}x{Q}_{n-1-k}x-\mu +\eta S_{n-1}x+D_1\frac{{\partial }^{2\alpha }{S}_{n-1}x}{\partial x^{2\alpha }}{t}^{k\alpha },\text{(25)}& E_{n}x,t=\sum _{n=1}^{\infty }\frac{\Gamma n-1\alpha +1}{\Gamma n\alpha +1}\frac{\beta }{N}\sum _{k=0}^{n-1}{S}_{k}x{I}_{n-1-k}x+S_{k}x{I}_{n-1-k}^{\ast }x+S_{k}x{I}_{n-1-k}x-{\alpha }_1+{\alpha }_2+\mu +\nu E_{n-1}x+D_2\frac{{\partial }^{2\alpha }{E}_{2}x}{\partial x^{2\alpha }}{t}^{k\alpha },\text{(26)}& I_{n}x,t=\sum _{n=1}^{\infty }\frac{\Gamma n-1\alpha +1}{\Gamma n\alpha +1}{\alpha }_1{E}_{n-1}x-{\gamma }_1+\mu +{\delta }_1+{\omega }_1{I}_{n-1}x+D_3\frac{{\partial }^{2\alpha }{I}_{n-1}x}{\partial x^{2\alpha }}{t}^{k\alpha },\text{(27)}& I_n^{\ast }x,t=\sum _{n=1}^{\infty }\frac{\Gamma n-1\alpha +1}{\Gamma n\alpha +1}{\alpha }_2{E}_{n-1}x-{\gamma }_2+\mu +{\delta }_2+{\omega }_2{I}_{n-1}^{\ast }x+D_4\frac{{\partial }^{2\alpha }{I}_{n-1}^{\ast }x}{\partial x^{2\alpha }}{t}^{k\alpha },\text{(28)}& Q_{n}x,t=\sum _{n=1}^{\infty }\frac{\Gamma n-1\alpha +1}{\Gamma n\alpha +1}\nu E_{n-1}x+{\omega }_1{I}_{n-1}x+{\omega }_2{I}_{n-1}^{\ast }x-\mu +\psi +{\delta }_3{Q}_{n-1}x{t}^{k\alpha },\text{(29)}& R_{n}x,t=\sum _{n=1}^{\infty }\frac{\Gamma n-1\alpha +1}{\Gamma n\alpha +1}\varphi Q_{n-1}x+{\gamma }_1{I}_{n-1}x+{\gamma }_2{I}_{n-1}^{\ast }x+\eta S_{n-1}x-\mu +\kappa R_{n-1}x{t}^{k\alpha }.\end{array}$$

3.2.1. Existence and Uniqueness of the Series Solution of the Nonsystem of Fractional PDEs

The proofs of the existence and uniqueness of the series solutions in equations (24)–(29) of the nonlinear system of fractional PDEs are provided therein. $$\begin{array}{c}\text{(30)}& P_{1}x,t,st=\frac{\Gamma n-1\alpha +1}{\Gamma n\alpha +1}\psi N+\kappa R_{n-1}x-\frac{\beta }{N}\sum _{k=0}^{n-1}{S}_{k}x{I}_{n-1-k}x+S_{k}x{I}_{n-1-k}^{\ast }x+S_{k}x{Q}_{n-1-k}x-\mu +\eta S_{n-1}x+D_1\frac{{\partial }^{2\alpha }{S}_{n-1}x}{\partial x^{2\alpha }}{t}^{k\alpha },\text{(31)}& P_{1}x,t,s^{\prime }t=\frac{\Gamma n-1\alpha +1}{\Gamma n\alpha +1}\psi N+\kappa R_{n-1}x-\frac{\beta }{N}\sum _{k=0}^{n-1}{S}_k^{\prime }x{I}_{n-1-k}x+S_k^{\prime }x{I}_{n-1-k}^{\ast }x+S_k^{\prime }x{Q}_{n-1-k}x-\mu +\eta S_{n-1}^{\prime }x+D_1\frac{{\partial }^{2\alpha }{S}_{n-1}^{\prime }x}{\partial x^{2\alpha }}{t}^{k\alpha },P_{1}x,t,st-Px,t,s^{\prime }=\frac{\Gamma n-1\alpha +1}{\Gamma n\alpha +1}-\frac{\beta }{N}\sum _{k=0}^{n-1}{I}_{n-1-k}x{S}_{k}x-S_k^{\prime }x+I_{n-1-k}^{\ast }x{S}_{k}x-S_k^{\prime }x+Q_{n-1-k}x{S}_{k}x-S_k^{\prime }x-\mu +\eta S_{n-1}x-S_{n-1}^{\prime }x+D_1\frac{{\partial }^{2\alpha }{S}_{n-1}x-S_{n-1}^{\prime }x}{\partial x^{2\alpha }}{t}^{k\alpha }\le \frac{\Gamma n-1\alpha }{\Gamma n\alpha +1}\frac{\beta }{N}\sum _{k=0}^{n-1}{I}_{n-1-k}x+I_{n-1-k}^{\ast }x+Q_{n-1-k}x{S}_{k}x-S_k^{\prime }x+\mu +\eta +D_1\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{S}_{n-1}x-S_{n-1}^{\prime }x{t}^{k\alpha }\le \frac{\Gamma n-1\alpha }{\Gamma n\alpha +1}\frac{\beta }{N}\sum _{k=0}^{n-1}{I}_{n-1-k}x+I_{n-1-k}^{\ast }x+Q_{n-1-k}x+\mu +\eta +D_1\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{S}_{n-1}x-S_{n-1}^{\prime }x{t}^{k\alpha },\\ P_{1}x,t,st-P_{1}x,t,s^{\prime }t\le \frac{\beta }{N}\frac{\Gamma n-1\alpha +1}{\Gamma n\alpha +1}{l}_1+l_2+l_3+\mu +\eta +D\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{S}_{k}x-S_k^{\prime }x,\\ P_{1}x,t,st-P_{1}x,t,s^{\prime }t\le {\lambda }_1{S}_{k}x-S_k^{\prime }x,\end{array}$$where $$\begin{array}{}\text{(32)}& {\lambda }_1=\frac{\beta }{N}\frac{\Gamma n-1\alpha +1}{\Gamma n\alpha +1}{l}_1+l_2+l_3+l_4+l_5,l_1=\sum _{k=0}^{n-1}{I}_{n-1-k}x,l_2=\sum _{k=0}^{n-1}{I}_{n-1-k}^{\ast }x,l_3=\sum _{k=0}^{n-1}{Q}_{n-1-k}x,l_4=\mu +\eta ,l_5=D_1\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }},\end{array}$$with $0<{\lambda }_1\le 1$.

This implies that the function is Lipschitz continuous on the domain $x,t,stx\in {\mathfrak{R}}^{+},t\in 0\cup {\mathfrak{R}}^{+}\text{and}st\in {\mathfrak{R}}^{+}$.

Following similar procedure above, the following continuous functions are obtained over the domain: $$\begin{array}{c}\text{(33)}& P_{2}x,t,st-P_{2}x,t,s^{\prime }t\le {\lambda }_2{E}_{k}x-E_k^{\prime },P_{2}x,t,Ex-P_{2}x,t,E^{\prime }x\le \frac{\Gamma n-1\alpha +1}{\Gamma n\alpha +1}\frac{\beta }{N}\sum _{k=0}^{n-1}{I}_{n-1-k}x+I_{n-1-k}^{\ast }x+Q_{n-1-k}x+{\alpha }_1+{\alpha }_2+\mu +\nu +D_2\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{E}_{n}x-E_n^{\prime }x,\\ Px,t,Ex-P_{2}x,t,E^{\prime }x\le {\lambda }_2{E}_{n}x-E_n^{\prime }x,\end{array}$$where $$\begin{array}{}\text{(34)}& {\lambda }_2=\frac{\beta }{N}\frac{\Gamma n-1\alpha +1}{\Gamma n\alpha +1}{l}_1+l_2+l_3+l_4+l_5+l_6,\end{array}$$where $l_6={\alpha }_1+{\alpha }_2+\mu +\nu$ and $l_1+l_2+l_3$ and $l_5$ have usual meanings. $$\begin{array}{}\text{(35)}& \begin{array}{c}Px,t,Ix-P_{3}x,t,I^{\prime }x\le {\lambda }_3{I}_{n-1}x-I_{n-1}^{\prime },\end{array}\end{array}$$where $$\begin{array}{}\text{(36)}& {\lambda }_3=\frac{\beta }{N}\frac{\Gamma n-1\alpha +1}{\Gamma n\alpha +1}{l}_1+l_2+l_3+l_4+l_5+l_6,\end{array}$$where, $l_6={\alpha }_1+{\alpha }_2+\mu +\nu$ and $l_1+l_2+l_3$ and $l_5$ have usual meanings.

3.3. Analytic Solutions of the System of Fractional Partial Differential Equations Using the Residual Power Series Method

This section contains the series solutions of the nonlinear system of equations (2)–(7), together with the initial conditions in equation (8). In using the RPSM, it is assumed that there are discrepancies between the terms on the right hand sides and the left hand sides of the system of equations (2)–(7). With this assumption, the approximations of the dependent variable with respect to only one independent variable are obtained depending on the given initial condition or boundary point condition. The other independent variable is automatically in fractional form which converges to a point in the Holder’s spaces. In doing this, we set $$\begin{array}{}\text{(37)}& S_{k}x,t=S_{o}x+\sum _{n=1}^k\frac{S_{n}x}{\Gamma n\alpha +1}{t}^{n\alpha },\text{(38)}& E_{k}x,t=E_{o}x+\sum _{n=1}^k\frac{E_{n}x}{\Gamma n\alpha +1}{t}^{n\alpha },\text{(39)}& I_{k}x,t=I_{o}x+\sum _{n=1}^k\frac{I_{n}x}{\Gamma n\alpha +1}{t}^{n\alpha },\text{(40)}& I_k^{\ast }x,t=I_o^{\ast }x+\sum _{n=1}^k\frac{I_n^{\ast }x}{\Gamma n\alpha +1}{t}^{n\alpha },\text{(41)}& Q_{k}x,t=Q_{o}x+\sum _{n=1}^k\frac{Q_{n}x}{\Gamma n\alpha +1}{t}^{n\alpha },\text{(42)}& R_k^{\ast }x,t=R_{o}x+\sum _{n=1}^k\frac{R_{n}x}{\Gamma n\alpha +1}{t}^{n\alpha }.\text{(43)}& \mathrm{Re}{s}_{k}x,t=\frac{{\partial }^{\alpha }S}{\partial t^{\alpha }}-\psi N-\kappa Rx,t+\frac{\beta }{N}Sx,tIx,t+I^{\ast }x,t+Qx,t+\mu +\eta Sx,t-D_1\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}Sx,t.\end{array}$$

Substituting equation (37), (39), (40), (41), and (42) into equation (43) yields $$\begin{array}{}\text{(44)}& \mathrm{Re}{s}_{k}x,t=\frac{{\partial }^{\alpha }}{\partial t^{\alpha }}{S}_{o}x+\sum _{n=1}^k\frac{S_{n}x}{\Gamma n\alpha +1}{t}^{n\alpha }-\psi N-\kappa R_{o}x+\sum _{n=1}^k\frac{R_{o}x}{\Gamma n\alpha +1}{t}^{n\alpha }+\frac{\beta }{N}{S}_{o}x+\sum _{n=1}^k\frac{S_{n}x}{\Gamma n\alpha +1}{t}^{n\alpha }{I}_{o}x+\sum _{n=1}^k\frac{I_{n}x}{\Gamma n\alpha +1}{t}^{n\alpha }+I_o^{\ast }x+\sum _{n=1}^k\frac{I_n^{\ast }x}{\Gamma n\alpha +1}{t}^{n\alpha }+Q_{o}x+\sum _{n=1}^k\frac{Q_{n}x}{\Gamma n\alpha +1}{t}^{n\alpha }+\mu +\eta S_{o}x+\sum _{n=1}^k\frac{S_{n}x}{\Gamma n\alpha +1}{t}^{n\alpha }-D_1\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{S}_{o}x+\sum _{n=1}^k\frac{S_{n}x}{\Gamma n\alpha +1}{t}^{n\alpha }.\end{array}$$

To obtain $S_{1}x$, equation (43) is reduced to $$\begin{array}{}\text{(45)}& \mathrm{Re}{s}_{1}x,t=\frac{{\partial }^{\alpha }}{\partial t^{\alpha }}{S}_{o}x+\frac{S_{1}x}{\Gamma \alpha +1}{t}^{\alpha }-\psi N-\kappa R_{o}x+\frac{R_{1}x}{\Gamma \alpha +1}{t}^{\alpha }+\frac{\beta }{N}{S}_{o}x+\frac{S_{1}x}{\Gamma \alpha +1}{t}^{\alpha }{I}_{o}x+\frac{I_{1}x}{\Gamma \alpha +1}{t}^{\alpha }+I_o^{\ast }x+\frac{I_1^{\ast }x}{\Gamma \alpha +1}{t}^{\alpha }+Q_{o}x+\frac{Q_{1}x}{\Gamma \alpha +1}{t}^{\alpha }+\mu +\eta S_{o}x+\frac{S_{1}x}{\Gamma \alpha +1}{t}^{\alpha }-D_1\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{S}_{o}x+\frac{S_{1}x}{\Gamma \alpha +1}{t}^{\alpha }.\end{array}$$

Setting $\mathrm{Re}sx,0=0$, it implies that $$\begin{array}{}\text{(46)}& S_{1}x=\psi N+\kappa R_{o}x-\frac{\beta }{N}{S}_{o}x{I}_{o}x+I_o^{\ast }x+Q_{o}x-\mu +\eta S_{o}x+D_1\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{S}_{o}x.\end{array}$$

Similarly, the $S_{2}x$ is obtained as $$\begin{array}{c}\text{(47)}& \mathrm{Re}{s}_{2}x,t=\frac{{\partial }^{\alpha }}{\partial t^{\alpha }}{S}_{0}x+\frac{S_{1}x}{\Gamma \alpha +1}{t}^{\alpha }+\frac{S_{2}x}{\Gamma 2\alpha +1}{t}^{2\alpha }-\psi N-\kappa R_{0}x+\frac{R_{1}x}{\Gamma \alpha +1}{t}^{\alpha }+\frac{R_{2}x}{\Gamma 2\alpha +1}{t}^{2\alpha }+\frac{\beta }{N}{S}_{0}x+\frac{S_{1}x}{\Gamma \alpha +1}{t}^{\alpha }+\frac{S_{2}x}{\Gamma 2\alpha +1}{t}^{2\alpha }{I}_{0}x+\frac{I_{1}x}{\Gamma \alpha +1}{t}^{\alpha }+\frac{I_{2}x}{\Gamma 2\alpha +1}{t}^{2\alpha }+I_0^{\ast }x+\frac{I_1^{\ast }x}{\Gamma \alpha +1}{t}^{\alpha }+\frac{I_2^{\ast }x}{\Gamma 2\alpha +1}{t}^{2\alpha }+Q_{0}x+\frac{Q_{1}x}{\Gamma \alpha +1}{t}^{\alpha }+\frac{Q_{2}x}{\Gamma 2\alpha +1}{t}^{2\alpha }+\mu +\eta S_{0}x+\frac{S_{1}x}{\Gamma \alpha +1}{t}^{\alpha }+\frac{S_{2}x}{\Gamma 2\alpha +1}{t}^{2\alpha }-D_1\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{S}_{0}x+\frac{S_{1}x}{\Gamma \alpha +1}{t}^{\alpha }+\frac{S_{2}x}{\Gamma 2\alpha +1}{t}^{2\alpha },\frac{{\partial }^{\alpha }}{\partial t^{\alpha }}\mathrm{Re}{s}_{2}x,0=S_{2}x-\kappa R_{1}x+\frac{\beta }{N}{S}_{1}x{I}_{0}x+I_0^{\ast }x+Q_{0}x+S_{0}x{I}_{1}x+I_1^{\ast }x+Q_{1}x+\mu +\eta S_{1}x-D_1\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{S}_{1}x=0,\\ S_{2}x=\kappa R_{1}x-\frac{\beta }{N}{S}_{1}x{I}_{0}x+I_0^{\ast }x+Q_{0}x+S_{0}x{I}_{1}x+I_1^{\ast }x+Q_{1}x-\mu +\eta S_{1}x+D_1\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{S}_{1}x.\end{array}$$

Similarly, the following results are obtained: $$\begin{array}{c}\text{(48)}& \frac{{\partial }^{\alpha }}{\partial t^{\alpha }}\mathrm{Re}{s}_{3}x,0=S_{3}x-\kappa R_{2}x+\frac{\beta }{N}{S}_{2}x{I}_{0}x+I_0^{\ast }x+Q_{0}x+2S_{1}x{I}_{1}x+I_1^{\ast }x+Q_{1}x+S_{0}x{I}_{2}x+I_2^{\ast }x+Q_{2}x+\mu +\eta S_{2}x-D_1\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{S}_{2}x=0,S_{3}x=\kappa R_{2}x-\frac{\beta }{N}{S}_{2}x{I}_{0}x+I_0^{\ast }x+Q_{0}x+2S_{1}x{I}_{1}x+I_1^{\ast }x+Q_{1}x+S_{0}x{I}_{2}x+I_2^{\ast }x+Q_{2}x-\mu +\eta S_{2}x+D_1\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{S}_{2}x\\ ⋮\\ S_{n}x=\kappa R_{n-1}x-\frac{\beta }{N}{S}_{n-1}x{I}_{0}x+I_0^{\ast }x+Q_{0}x+\begin{array}{c}n-1\\ 1\end{array}{S}_{n-2}x{I}_{1}x+I_1^{\ast }x+Q_{1}x+\cdots +\begin{array}{c}n-1\\ r\end{array}{S}_{n-1-r}x{I}_{r}x+I_r^{\ast }x+Q_{r}x+\cdots +S_{0}x{I}_{n-1}x+I_{n-1}^{\ast }x+Q_{n-1}x-\mu +\eta S_{n-1}x+D_1\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{S}_{n-1}x,\\ S_{n}x=\kappa R_{n-1}x-\frac{\beta }{N}\sum _{r=0}^{n-1}\begin{array}{c}n-1\\ r\end{array}{S}_{n-1-r}x{I}_{r}x+I_r^{\ast }x+Q_{r}x-\mu +\eta S_{n-1}x+D_1\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{S}_{n-1}x.\end{array}$$

Similarly, we obtain the following results for $Ex$ as $$\begin{array}{c}\text{(49)}& E_{1}x=\frac{\beta }{N}{S}_{0}x{I}_{0}x+I_0^{\ast }x+Q_{0}x-{\alpha }_1+{\alpha }_2+\mu +\nu E_{0}x+D_2\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{E}_{0}x,E_{2}x=\frac{\beta }{N}{S}_{1}x{I}_{0}x+I_0^{\ast }x+Q_{0}x+S_{0}x{I}_{1}x+I_1^{\ast }x+Q_{1}x-{\alpha }_1+{\alpha }_2+\mu +\nu E_{1}x+D_2\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{E}_{1}x,\\ E_{3}x=\frac{\beta }{N}{S}_{2}x{I}_{0}x+I_0^{\ast }x+Q_{0}x+2S_{1}x{I}_{1}x+I_1^{\ast }x+Q_{1}x+S_{0}x{I}_{2}x+I_2^{\ast }x+Q_{2}x-{\alpha }_1+{\alpha }_2+\mu +\nu E_{2}x+D_2\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{E}_{2}x,\\ ⋮\\ E_{n}x=\frac{\beta }{N}\sum _{r=0}^{n-1}\begin{array}{c}n-1\\ r\end{array}{S}_{n-1-r}x{I}_{r}x+I_r^{\ast }x+Q_{r}x-{\alpha }_1+{\alpha }_2+\mu +\nu E_{n-1}x+D_2\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{E}_{n-1}x.\end{array}$$

Following similar procedure, the series solutions of the fractional PDEs for describing the number of deviant infectives, number of non-deviant infectives, number of quarantined persons, and number of recoveries from the COVID-19 disease are obtained as follows: $$\begin{array}{c}\text{(50)}& I_{1}x={\alpha }_1{E}_{0}x-{\gamma }_1+\mu +{\delta }_1+{\omega }_1{I}_{0}x+D_3\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{I}_{0}x,I_{2}x={\alpha }_1{E}_{1}x-{\gamma }_1+\mu +{\delta }_1+{\omega }_1{I}_{1}x+D_3\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{I}_{1}x,\\ I_{3}x={\alpha }_1{E}_{2}x-{\gamma }_1+\mu +{\delta }_1+{\omega }_1{I}_{2}x+D_3\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{I}_{2}x,\\ ⋮\\ I_{n}x={\alpha }_1{E}_{n-1}x-{\gamma }_1+\mu +{\delta }_1+{\omega }_1{I}_{n-1}x+D_3\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{I}_{n-1}x,\\ I_1^{\ast }x={\alpha }_2{E}_{0}x-{\gamma }_2+\mu +{\delta }_2+{\omega }_2{I}_0^{\ast }x+D_4\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{I}_0^{\ast }x,\\ I_2^{\ast }x={\alpha }_2{E}_{1}x-{\gamma }_2+\mu +{\delta }_2+{\omega }_2{I}_1^{\ast }x+D_4\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{I}_1^{\ast }x,\\ I_3^{\ast }x={\alpha }_2{E}_{2}x-{\gamma }_2+\mu +{\delta }_2+{\omega }_2{I}_2^{\ast }x+D_4\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{I}_2^{\ast }x,\\ ⋮\\ I_n^{\ast }x={\alpha }_2{E}_{n-1}x-{\gamma }_2+\mu +{\delta }_2+{\omega }_2{I}_{n-1}^{\ast }x+D_4\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{I}_{n-1}^{\ast }x,\\ Q_{1}x=v{E}_{0}x+{\omega }_1{I}_{0}x+{\omega }_2{I}_0^{\ast }x-\mu +\varphi +{\delta }_3{I}_{0}x,\\ Q_{2}x=v{E}_{1}x+{\omega }_1{I}_{1}x+{\omega }_2{I}_1^{\ast }x-\mu +\varphi +{\delta }_3{I}_{1}x,\\ Q_{3}x=v{E}_{2}x+{\omega }_1{I}_{2}x+{\omega }_2{I}_2^{\ast }x-\mu +\varphi +{\delta }_3{I}_{2}x,\\ ⋮\\ Q_{n}x=v{E}_{n-1}x+{\omega }_1{I}_{n-1}x+{\omega }_2{I}_{n-1}^{\ast }x-\mu +\varphi +{\delta }_3{I}_{n-1}x,\\ R_{1}x=\varphi Q_{0}x+{\gamma }_1{I}_{0}x+{\gamma }_2{I}_0^{\ast }x+\eta S_{0}x-\mu +\kappa R_{0}x,\\ R_{2}x=\varphi Q_{1}x+{\gamma }_1{I}_{1}x+{\gamma }_2{I}_1^{\ast }x+\eta S_{1}x-\mu +\kappa R_{1}x,\\ R_{3}x=\varphi Q_{2}x+{\gamma }_1{I}_{2}x+{\gamma }_2{I}_2^{\ast }x+\eta S_{2}x-\mu +\kappa R_{2}x,\\ ⋮\\ R_{n}x=\varphi Q_{n-1}x+{\gamma }_1{I}_{n-1}x+{\gamma }_2{I}_{n-1}^{\ast }x+\eta S_{n-1}x-\mu +\kappa R_{n-1}x.\end{array}$$

In order to obtain the series solution for the number of susceptible individuals, we substitute the last equation of the system of equations (48) into equation (37) which yields $$\begin{array}{}\text{(51)}& S_{k}x,t=S_{0}x+\sum _{n=1}^k\frac{1}{\Gamma n\alpha +1}\kappa R_{n-1}x-\frac{\beta }{N}\sum _{r=0}^{n-1}\begin{array}{c}n-1\\ r\end{array}{S}_{n-1-r}x{I}_{r}x+I_r^{\ast }x+Q_{r}x-\mu +\eta S_{n-1}x+D_1\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{S}_{n-1}x{t}^{n\alpha }.\end{array}$$

Similarly, the following results are obtained for the number of exposed persons, deviant infectives, nondeviant infectives, quarantined persons, and recovered persons: $$\begin{array}{c}\text{(52)}& E_{k}x,t=E_{0}x+\sum _{n=1}^k\frac{1}{\Gamma n\alpha +1}\frac{\beta }{N}\sum _{r=0}^{n-1}\begin{array}{c}n-1\\ r\end{array}{S}_{n-1-r}x{I}_{r}x+I_r^{\ast }x+Q_{r}x-{\alpha }_1+{\alpha }_2+\mu +\nu E_{n-1}x+D_2\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{E}_{n-1}x{t}^{n\alpha },I_{k}x,t=I_{0}x+\sum _{n=1}^k\frac{1}{\Gamma n\alpha +1}{\alpha }_1{E}_{n-1}x-{\gamma }_1+\mu +{\delta }_1+{\omega }_1{I}_{n-1}x+D_3\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{I}_{n-1}x{t}^{n\alpha },\\ \text{(53)}& I_k^{\ast }x,t=I_0^{\ast }x+\sum _{n=1}^k\frac{1}{\Gamma n\alpha +1}{\alpha }_2{E}_{n-1}x-{\gamma }_2+\mu +{\delta }_2+{\omega }_2{I}_{n-1}^{\ast }x+D_4\frac{{\partial }^{2\alpha }}{\partial x^{2\alpha }}{I}_{n-1}^{\ast }x{t}^{n\alpha },\text{(54)}& Q_{k}x,t=Q_{0}x+\sum _{n=1}^k\frac{1}{\Gamma n\alpha +1}v{E}_{n-1}x+{\omega }_1{I}_{n-1}x+{\omega }_2{I}_{n-1}^{\ast }x-\mu +\varphi +{\delta }_3{I}_{n-1}x{t}^{n\alpha },\text{(55)}& R_{k}x,t=R_{0}x+\sum _{n=1}^k\frac{1}{\Gamma n\alpha +1}\varphi Q_{n-1}x+{\gamma }_1{I}_{n-1}x+{\gamma }_2{I}_{n-1}^{\ast }x+\eta S_{n-1}x-\mu +\kappa R_{n-1}x{t}^{n\alpha }.\end{array}$$

3.4. Numerical Results

In this section, the three-dimensional plots, as well as the two-dimensional plots, are provided here. The initial condition for each subgroup was taken from [28].

3.4.1. The Numerical Results for the Series Solutions of the Nonlinear System of Fractional PDEs Using the FPSM

This section contains the plots for the series solutions to the system of equations (22)–(27) for $\alpha =0.2000$. Both three-dimensional and two-dimensional plots for each series solution of these equations are depicted in Figures 2 and 3, respectively. The initial conditions for each subgroup of the population were estimated using the available data in [28]. The number of susceptible decreases rapidly as they come into contact with the infectives, both locals and foreigners, as shown in Figure 2(a). In contrast to the number of nondeviant infective persons (see Figure 3(a)), the number of deviant infectives unexpectedly grows to its peak for a long time before it starts to decline and eventually becomes asymptotically stable on the $t$-axis. This implies that patients with the intention of transmitting COVID-19 disease to the susceptible members are the primary disease spreaders in Ghana.

[figure(s) omitted; refer to PDF]

For the number of susceptible, exposed, deviant, nondeviant, quarantined, and recovered subgroups, spatial distance is suppressed and time is varied. The first three series solutions for the system of equations (51)-(55) were employed to show trends in the numbers of susceptible individuals, exposed individuals, deviant infectives, nondeviant infectives, quarantined individuals, and recovered individuals. The plots were repeated for the first fifteen terms of the series solutions of system of equations (51)-(55), shown in Figure 4.

[figure(s) omitted; refer to PDF]

To take into account the impact of each series as the number of terms rises, plots for the first 80 terms of the series solutions of the system of equations (51)-(55) were reproduced as displayed in Figure 5. While the curves for the quarantined and recovered populations climb rapidly, the curve for nondeviant infected individuals declines sharply. This indicates that the majority of nondeviant carriers of the COVID-19 virus are not infecting susceptible members who are at risk of becoming infected. On the other hand, a large number of COVID-19 patients are being isolated at numerous facilities across the country, and these people are making progress every day. The fact that the deviant subpopulation is still growing suggests that patients have been spreading the SARS virus for a considerable amount of time with the intention of infecting a sizable number of vulnerable individuals with the COVID-19 infection. The plot for the exposed individuals increases to a peak and then declines and asymptotically moves toward the $t$-axis, showing that those who are vulnerable to contracting COVID-19 disease are at high risk. The SARS virus starts spreading to anybody who come into contact with it after a brief time of incubation. From day zero, the susceptible curve declines and becomes asymptotically toward the $t$-axis, signalling the end of the disease.

[figure(s) omitted; refer to PDF]

In Figure 6, there is relatively positive nonlinear relationship between the susceptible subgroups. This indicates that susceptible (foreigners) continuous inflow into the country irrespective of their home country will not ensue to the epidemiology of COVID-19 disease in Ghana. However, there is a slight positive nonlinear relationship and fairly negative nonlinear relationship between the nondeviant subpopulation and the distance.

[figure(s) omitted; refer to PDF]

3.4.2. The Numerical Results for the Series Solutions of the Nonlinear System of Fractional PDEs Using the RPSM

Figure 7 displays the first three terms of the series solutions, equations (51)-(55), with $\alpha =0.2000$. Every subgroup of the population size in Ghana corresponds with the epidemiological pattern of the COVID-19 disease.

[figure(s) omitted; refer to PDF]

3.5. Comparison of the FPSM and the RPSM

This section uses quantitative results of the FPSM to compare the quantitative results of the RPSM which are then superimposed on the quantitative results from the field data.

Figure 8 displays the series solutions of the first three terms using both the FPSM and the RPSM. The FPSM shows that the solution rises from the starting point to the peak and then falls, showing the loss of the susceptible members to the exposed subgroup and the exposed subgroup’s loss of members to both the deviant and nondeviant subgroups. On the other hand, the length of the series solution obtained using the RPSM increases starting with the beginning of the COVID-19 pandemic and takes a lot of time. It begins to decline in the direction of the $t$-axis, showing that the susceptible individuals who contract the SARS virus also persistently infect the population subjects. Despite this, the series solution using the FPSM is more consistent with the field data at the start of the disease outbreak than the series solution using the RPSM.

[figure(s) omitted; refer to PDF]

On the other hand, using the RPSM is more consistent with the field data as the COVID-19 infection is present in the population subjects over a long period of time in comparison to the FPSM’s series solution. This is indicated by the deep red plot. A similar observation was made in Figure 9 when comparing the two series solutions using both the FPSM and the RPSM.

[figure(s) omitted; refer to PDF]

4. Discussion

In contrast, the FPSM series solution for the number of susceptible individuals is proportional to the RPSM series solution for the number of susceptible individuals, with a proportional constant of $\psi \Gamma n-1\alpha +1$. The series solution of the vulnerable members utilizing the RPSM then reduced quickly in comparison to the series solution using the FPSM. For instance, the third term of the series solution using the RPSM is reduced by $S_{1}x{I}_{1}x+I_1^{\ast }x+Q_{1}x$ when compared to the series solution using the FPSM. In comparison, using the RPSM in obtaining the series solutions of the nonlinear system of equations, (2)–(7), together with the initial conditions in equation (8), are more consistent with the field data as compared to the series solutions given by the FPSM, as displayed in Figures 8 and 9.

5. Conclusion

There are more terms for the total number of susceptible members and exposed individuals in the series of solutions of the nonlinear system of fractional PDEs provided by the RPSM than the series of solutions yielded by the FPSM. The nonlinear term $SxIx+I^{\ast }x+Qx$ is what causes the difference between the two series solutions for the number of susceptible members and exposed individuals. The power of $SxIx+I^{\ast }x+Qx$ is linear for the first and second terms and then increases in the pattern of natural numbers, that is, having positive integer binomial powers, when employing the RPSM. Because the product of this nonlinear term is differentiated, this occurs. The nonlinear term, however, increases linearly when the FPSM is used to find the series solutions of the system of equations (2)–(7) along with the initial conditions in equation (8). Interestingly, series solutions of the nonlinear system of the equations using the RPSM were observed to be more consistent as compared to series solutions given by FPSM. This is due to the fact that the RPSM utilizes the variations in the nonlinear system of equations unlike the FPSM. However, both the RPSM and the FPSM yield the same series solutions of the linear system of equations, as indicated by equations (24)–(29) and system of equations (51)–(55).

References

[1] I. A. Baba, B. A. Baba, P. Esmaili, "A mathematical model to study the effectiveness of some of the strategies adopted in curtailing the spread of COVID-19," Computational and Mathematical Methods in Medicine, vol. 2020,DOI: 10.1155/2020/5248569, 2020.

[2] J. Y. Mugisha, J. Ssebuliba, J. N. Nakakawa, C. R. Kikawa, A. Ssematimba, "Mathematical modeling of COVID-19 transmission dynamics in Uganda: implications of complacency and early easing of lockdown," PLoS One, vol. 16 no. 2, article e0247456,DOI: 10.1371/journal.pone.0247456, 2021.

[3] B. Barnes, J. Ackora-Prah, F. O. Boateng, L. Amanor, "Mathematical modelling of the epidemiology of COVID-19 infection in Ghana," Scientific African, vol. 15, article e01070,DOI: 10.1016/j.sciaf.2021.e01070, 2022.

[4] WHO, "WHO director-general’s opening remarks at the media briefing on COVID-19. World Health Organization," 2020. http://www.who.int/director-general/speeches/detail/who-director-general-s-opening-remarks-at-the-media-briefing-on-covid-19-11-match-2020

[5] A. Viguerie, A. Veneziani, G. Lorenzo, D. Baroli, N. Aretz-Nellesen, A. Patton, T. E. Yankeelov, A. Reali, T. J. R. Hughes, F. Auricchio, "Diffusion–reaction compartmental models formulated in a continuum mechanics framework: application to COVID-19, mathematical analysis, and numerical study," Computational Mechanics, vol. 66 no. 5, pp. 1131-1152, DOI: 10.1007/s00466-020-01888-0, 2020.

[6] B. Barnes, I. Takyi, B. E. Owusu, F. O. Boateng, A. Saahene, E. S. Baidoo, J. A. Adombire, "Mathematical modelling of the spatial epidemiology of COVID-19 with different diffusion coefficients," Journal of Differential Equations, vol. 2022, article 7563111,DOI: 10.1155/2022/7563111, 2022.

[7] J. O. Watson, G. Barnsley, J. Toor, B. A. Hogan, P. Winskill, C. A. Ghani, "Global impact of the first year of COVID- 19 vaccination: a mathematical modelling study," The Lancet Infectious Diseases, vol. 22 no. 9, pp. 1293-1302, DOI: 10.1016/S1473-3099(22)00320-6, 2022.

[8] A. M. Elaiw, R. S. Alsulami, A. D. Hobiny, "Modeling and stability analysis of within-host IAV/SARS-CoV-2 coinfection with antibody immunity," Mathematics, vol. 10 no. 22,DOI: 10.3390/math10224382, 2022.

[9] A. D. Algarni, A. B. Hamed, M. Hamed, M. Hamdi, H. Elmannai, S. Meshoul, "Mathematical COVID-19 model with vaccination. A case study in Saudi Arabia," PeerJ Computer Science, vol. 8, article e959,DOI: 10.7717/peerj-cs.959, 2022.

[10] M. L. Diagne, H. Rwezaura, S. Y. Tchoumi, J. M. Tchuenche, "A mathematical model of COVID-19 with vaccination and treatment," Computational and Mathematical Methods in Medicine, vol. 2021,DOI: 10.1155/2021/1250129, 2021.

[11] M. Udriste, M. Ferrara, D. Zugravescu, F. Munteanu, "Controllability of a nonholonomic macroeconomic system," Journal of Optimization Theory and Applications, vol. 154 no. 3, pp. 1036-1054, DOI: 10.1007/s10957-012-0021-x, 2012.

[12] I. Mohammed, A. Nauman, P. Paul, S. Ganesan, K. H. Chen, S. M. S. Jalil, D. Zakaria, "The efficacy and effectiveness of the COVID-19 vaccines in reducing infection, severity, hospitalization, and mortality: a systematic review," Human Vaccines and Immunotherapeutics, vol. 18 no. 1, article 2027160,DOI: 10.1080/21645515.2022.2027160, 2022.

[13] M. Shah, M. Arfan, I. Mahariq, A. Ahmadian, S. Salahshour, M. Ferrara, "Fractal-fractional mathematical model addressing the situation of corona virus in Pakistan," Results in Physics, vol. 19, article 103560,DOI: 10.1016/j.rinp.2020.103560, 2020.

[14] H. Vazquez-Leal, B. Benhammouda, U. Filobello-Nino, A. Sarmiento-Reyes, V. M. Jimenez-Fernandez, J. L. Garcia-Gervacio, J. Huerta-Chua, L. J. Morales-Mendoza, M. Gonzalez-Lee, "Direct application of Padé approximant for solving nonlinear differential equations," Springerplus, vol. 3 no. 1,DOI: 10.1186/2193-1801-3-563, 2014.

[15] O. A. Arqub, "Series solution of fuzzy differential equations under strongly generalized differentiability," Journal of Advanced Research in Applied Mathematics, vol. 5 no. 1, pp. 31-52, DOI: 10.5373/jaram.1447.051912, 2013.

[16] M. Modanli, S. T. Abdulazeez, A. M. Husien, "A residual power series method for solving pseudo hyperbolic partial differential equations with nonlocal conditions," Numerical Methods for Partial Differential Equations, vol. 37 no. 3, pp. 2235-2243, DOI: 10.1002/num.22683, 2021.

[17] A. El-Ajou, O. A. Arqub, S. Momani, D. Baleanu, A. Alsaedi, "A novel expansion iterative method for solving linear partial differential equations of fractional order," Applied Mathematics and Computation, vol. 257, pp. 119-133, DOI: 10.1016/j.amc.2014.12.121, 2015.

[18] M. A. Bayrak, A. Demir, E. Ozbilge, "An improved version of residual power series method for space-time fractional problems," Advances in Mathematical Physics, vol. 2022,DOI: 10.1155/2022/6174688, 2022.

[19] B. Chen, L. Qin, F. Xu, J. Zu, "Applications of general residual power series method to differential equations with variable coefficients," Discrete Dynamics in Nature and Society, vol. 2018,DOI: 10.1155/2018/2394735, 2018.

[20] P. Dunnimit, A. Wiwatwanich, D. Poltem, "Analytical solution of nonlinear fractional Volterra population growth model using the modified residual power series method," Symmetry, vol. 12 no. 11,DOI: 10.3390/sym12111779, 2020.

[21] S. Salahshour, A. Ahmadian, M. Salimi, M. Ferrara, D. Baleanu, "Asymptotic solutions of fractional interval differential equations with nonsingular kernel derivative," Chaos, vol. 29 no. 8, article 083110,DOI: 10.1063/1.5096022, 2019.

[22] S. Z. Rida, A. A. M. Arafa, "New method for solving linear fractional differential equations," International Journal of Differential Equations, vol. 2011,DOI: 10.1155/2011/814132, 2011.

[23] A. I. Ali, M. Kalim, A. Khan, "Solution of fractional partial differential equations using fractional power series method," International Journal of Differential Equations, vol. 2021,DOI: 10.1155/2021/6385799, 2021.

[24] M. Bagheri, A. Khani, "Analytical method for solving the fractional order generalized kdv equation by a beta-fractional derivative," Advances in Mathematical Physics, vol. 2020,DOI: 10.1155/2020/8819183, 2020.

[25] A. A. Kilbas, H. M. Srivastava, J. J. Trujillo, "Preface," Theory and Applications of Fractional Differential Equations,DOI: 10.1016/S0304-0208(06)80001-0, 2006.

[26] S. Das, Functional Fractional Calculus, vol. 1,DOI: 10.1007/978-3-642-20545-3, 2011.

[27] J. M. Morel, F. Takens, B. Teissier, The Analysis of Fractional Differential Equations, 2004.

[28] H. Ritchie, E. Mathieu, L. Rodés-Guirao, C. Appel, C. Giattino, E. Ortiz-Ospina, J. Hasell, B. Macdonald, D. Beltekian, M. Roser, "Coronavirus pandemic (COVID-19)," , 2020. https://ourworldindata.org/coronavirus