Application of Homotopy Analysis Method for Solving an SEIRS Epidemic Model

In this paper, we modified the model of [23] and then applied a new semi-analytic technique namely the Homotopy Analysis Method (HAM) in solving the SEIRS Epidemic Mathematical Model. The modified SEIRS model wasfirst formulated and adequately analyzed. We investigated the basic properties of the model by proving the positivity of the solutions and establishing the invariant region. We further obtained the steady states: disease-free equilibrium (DFE) and endemic equilibrium (EE), then we went further to determine the local stability of the DEF and EE using the basic reproduction number which was calculated. We also applied Lyaponuv method to prove the global stability of endemic equilibrium, The HAM was applied to obtain an accurate solution to the model in few iterations. Finally, a numerical solution (simulation) of the model was obtained using MAPLE 15 computation software.


Introduction
Mathematical models have been used in comparing, planning, implementing, evaluating and optimizing various detection, prevention, therapy and control programs.
World Health Statistics (2000) shown that some vector borne diseases as Malaria, Dengue and Yellow fever, continue to threaten throughout the public health of many people. Most of the biological problems are inherently nonlinear.
Modeling of natural phenomena in science and engineering mostly leads to nonlinear problems. Solutions of these nonlinear problems are more difficult compared to linear ones, especially through the analytical approach. Therefore, endless efforts are devoted either to find ways to solve them or to decrease the error in the solutions.
At present, there exist few well known analytical approaches for solving nonlinear problems, viz. perturbation technique; non-perturbation techniques (for example, Lyapunov's small parameter method, δ-expansion method, Adomian's decomposition method, etc.). Perturbation technique is widely applied in science and engineering because in this technique the original nonlinear equations are replaced by an infinite number of linear sub-problems which are easy to solve. Unfortunately, not every nonlinear even if there exists a small parameter, the sub-problem mighthave no solutions or might be rather complicated so that only a few of the sub-problems can be solved. Thus, it is not guaranteed that one can always apply perturbation approximations efficiently for a given nonlinear problem. Many new techniques have been proposed recently to eliminate the ''small parameter'' assumption, such as the artificial parameter method in [14], Variational iteration method in [6] and Homotopy analysis method in [13]. Among these methods, Homotopy Analysis Method (HAM) gives us an acceptable and more accurate analytical result with convenient convergence and stability. Homotopy Analysis Method and Modified Homotopy Perturbation Methodwaspresentedfor Strongly Nonlinear Oscillator. In his work he compared the analytical result as well as the numerical result obtained by using Homotopy Analysis Method (HAM), Homotopy Perturbation Method (HPM) and Modified Homotopy Perturbation Method (mHPM). In his study he used different examples with strongly nonlinearity to assess in detail and to illustrate the effectiveness and convenience of the methods. Comparing the three methods with attention focused on the accuracy of the results and applicability of the methods. He came to a conclusion that the solution using HAM which is valid in a small region is more accurate than the solution obtained using other two methods, [18].
Analytic technique, namely Homotopy Analysis Method (HAM) was used for solving an SEIR model. In his work HAM introduced a non-zero auxiliary parameter h to construct a two-parameter family of equations (the zerothorder deformation equation). The result of the theoretical analysis of the HAM in his work shows that it yields more accurate results in little iteration, [15].
A Homotopy Analysis Method was presented for SEIR tuberculosis model. In their paper, they provided a very accurate, non-perturbative, semi-analytical solution to a system of nonlinear first-order differential equations modeling the transmission of tuberculosis (TB) in a homogeneous population. Their analysis which is based on Homotopy Analysis Method (HAM) used Maple 15 Software to carry out their computations, and their result showed the validity and potential of HAM for computing the solution of nonlinear equations, [7].
An SEIR Model of Seasonal Epidemic Diseases using HAM was presented. The model was centered on childhood diseases measles, chickenpox mumps, rubella incorporate seasonal variation in contact rates due to the increased mixing during school terms compared to school holidays. Driven by seasonality the diseases were characterized by annual oscillations with variation contact rate which is a periodic function of time in years. They used HAM to obtain a semi analytic approximate solution of nonlinear simultaneous differential equations. Mathematical was used to carry out their computations. Result established through graphs showed the validity and potential of HAM for amplitude of Variation greater than zero, [17].
Homotopy analysis method (HAM) was applied to solve Non-linear Cauchy problem. They showed by of means HAM, that the solutions of some linear Cauchy problems of parabolic-hyperbolic type are exactly obtained in the form of convergent Taylor series. The HAM contains the auxiliary parameter, h, that provide a convenient way of controlling the convergent region of series solutions. The results reveal that the proposed method is very effective, [5].
Homotopy analysis method was employed in the analysis of one dimensional heat conduction problem with temperature dependent thermal conductivity. The temperature distribution is obtained in terms of the scale space variable x and as a function of B which is a parameter in the heat conduction equation. It was observed that the parameter B has a strong influence over the rate of heat conduction in the bar. By choosing the convergence parameter, h, in a suitable way, they obtained solutions for the temperature distribution for various values of B. From this temperature distribution other heat transfer quantities can be obtained. It is observed that for some values of B the temperature distribution in the bar is an increasing function of B while for other values of B it is a decreasing function. In their work, they introduced the basic idea of HAM and applied it to find an approximate analytical solution of the nonlinear one-dimensional heat conduction problem with temperature dependent thermal conductivity, [3].
Homotopy analysis method HAM was applied to solve Convection-Diffusion equation. In their work, the homotopy analysis method was used to find series solution of linear convection diffusion (CD) equation. They further explained that convection diffusion equation is a combination of the diffusion and convection (advection) equations and describes physical phenomena where particles, energy or other physical quantities are transferred inside a physical system due to two processes: diffusion and convection. This equation is a parabolic differential equation which describes physical phenomena. In their work they considered the following linear convection diffusion (CD) equation and the homotopy analysis method was applied to solve it, [4].
Hence, in their work they prove the validity and reliability of the homotopy analysis method (HAM) and with examples showed the accuracy and efficiency of the method. Therefore the homotopy analysis method is able to evaluate a satisfactory solution for convection diffusion equation (CD) Equations.
Homotopy Analysis Method was used for Non-linear periodic Oscillating equations with Absolute Value Term. Based on the homotopy analysis method (HAM), an analytic approach is proposed for highly nonlinear periodic oscillating equations with absolute value terms. The non-smoothness of absolute value terms is handled by means of Fourier expansion and the convergence is accelerated by using the iteration method. Two typical examples which cannot be solved by the method of averaging of perturbation technique are employed to illustrate the validity and flexibility of this approach. Rather accurate approximations are obtained using the HAM based approach. Hence the HAM approach has general meanings and thus can be used to solve many Absolute Valve term. This work coupled with Ramani equations were solved using Homotopy Analysis Method (HAM). The homotopy analysis method was successfully used to find out the approximate solutions of coupled Ramanni equations. The proposed method prepares one with an expedient way to control the convergence of approximate solution series. The basic ideas of this approach can be widely utilized to solve other strongly nonlinear evolution problems. Their numerical results indicated that the obtained approximate solutions were in suitable agreement with the exact solutions, thus demonstrating the remarkable efficiency of the homotopy analysis method.
Homotopy analysis method was applied for solving integral and Integro-differential equations. The homotopy analysis method was successfully applied to find the solution of integral and integro-differential equations. It is apparently seen that HAM is a very powerful and efficient technique in finding analytical solutions for wide classes of and integral equations. It is worth pointing out that this method presents a rapid convergence for the solutions. They also do not require large computer memory and discretization of the variables x. Their results show that HAM is a powerful mathematical tool for solving linear equations, [22].
Homotopy analysis method was applied for solving nonlinear Cauchy problem of parabolic-hyperbolic type are exactly obtained in the form of convergent Taylor series. The HAM contains the auxiliary parameter h that provides a convenient way of controlling the convergent region of series solutions. Their results revealed that the proposed method is very effective and simple, [5].
A research carried out that deals with mathematical modeling of malaria transmission in North Senatorial Zone of Taraba State, Nigeria. The SIR proposed by Kermack and McKendrick was implemented using data obtained from Essential Programme on Immunisation (EPI) unit, F.M.C., Jalingo, Taraba state, the data was used to analyse the rate of infection of malaria in the zone. From their analysis, they found out that the reproduction ratio is greater than 0, which implies that the force of malaria infection in Taraba North Senatorial Zone is high. It was recommended specific steps for the reduction of malaria in the zone, [15].
A researcher carried out a review study on mathematical models of malaria with emphasis on their relevance for control. They concluded that; Ross-Macdonald model of malaria transmission has had major influence on malaria control, and that one of its main contributions is that endemicity of malaria is most sensitive to changes in mosquito imago survival rate, and hence malaria can be controlled more efficiently with imagicites than with larvicides. It was further opined that Models that describe the immune response and simulate vaccination programs suggest that one of the most important determinants of the outcome of a vaccine campaign is the duration of vaccine efficacy, [11]. Apparently malaria can be controlled only if the duration of efficacy is in the order of a human life-span.
A mathematical model to study the transmission of malariawas developed; the authors showed from the study that the model has unique disease-free and endemic equilibria in the application of treatment and preventive measures. We analyze the existence and stability of disease-free and endemic malaria equilibria. Key to their study is the definition of a reproductive number R0 which is the number of the new infections caused by one individual in an otherwise fully susceptible population, through the duration of the infectious period. The disease-free equilibrium is locally asymptotically stable, if R0≤1 and that the endemic equilibrium exist provided R0>1. In fact, their study showed clearly that effective treatment offered to about fifty percent of the infected population together with about fifty percent prevention rate is all that is required to eliminate the disease, [1].
A model was developed and in their model they showed that the three states which a human can be in are; susceptible, infectious and recovered. However, in reality the recovered class is not removed or quarantined from the main population, in fact they become part of the susceptible class. In this work, they present a reviewed dynamical model for the malaria parasite/host dynamics on this basis, [10].
Other papers consulted and referenced are [2,8,9,12,16,19,20] In this paper, we modified the model of [23]and then solve it with HAM. The total population size N(t) is divided into four compartment, a susceptible compartment, labeled S,( in which all the individuals are vulnerable to the disease but are not yet infected with the disease), an exposed compartment, labeled E,( in which all the individuals are already infected with the disease but are not yet infectious, that is, they cannot infect others with the disease), an infected compartment, labeled I (in which all the individuals are infected by the disease and infectious), a recovered compartment, labeled R in which all the individuals are recovered from the disease and we assume that immunity is not permanent and that recovered individuals do revert to the susceptible class), hence joins the susceptible class again, S. Denote the numbers of individuals in the compartments S, E, I, R, and S, at time t as S(t), E(t), I(t), R(t), and S(t) respectively.

Flow Diagram
We first draw the flow diagram for [22]. We now present the flow diagram of the modified model.

Assumptions of the Model
1. For this model we assume that Immunity is not permanent i.e. Recovered individuals can become susceptible again. 2. We also assume that individuals enter the population by birth or immigration and exist by death or infection. 3. The incidence term is of the standard bilinear mass action this is due to the very large population of infected people and lastly we assume that individuals enter the population by birth or immigration with a constant recruitment as susceptible S and exit by death or by infection.

Symbols and Parameters
Here we show the description of the symbols and parameter used in formulating the model.

S(t) The number of susceptible individuals at time t E(t)
The number of exposed individuals at time t I(t) The number of Infectious individuals at time t R(t) The number of Recovered individuals at time t Λ (t) Recruitment rate β (t) Effective contact rate ε (t) The rate at which exposed individuals become infectious γ (t) The rate at which infected individuals recover ω (t) The rate at which recovered individuals become

Variables and Parameters Description
The total population number

Mathematical Model Formulation
Formulation of the SEIRS mathematical model using the flow diagram above we have:

Basic Properties of the Model
Note that In the absence of the disease, 0 In the absence of disease, the total population size N(t) converges to the equilibrium µ Λ .
From epidemiology consideration, we study our model (1) in the closed feasible region:

Positivity of Solutions
We show that all the variables in the model equation are non-negative.
Considering the second equation of (2.1) we have, Following the same way, we can prove that ( ) 0 ( ) 0 I t and R t ≥ ≥ Hence, for t> 0 S(t), E(t), I(t) and R(t) are all non-negative.

Existence of Steady State of the Model
At steady state This implies that Solving (2) we have The disease-free equilibrium is

The Model Basic Reproduction Number
The local stability is established by using the next generation operator method on the system (1).
The Basic Reproduction Number ℜ0 is defined as the expected number of secondary cases produced, in a completely susceptible population, by a typical infective individual [21]. This definition is given for the models that represent spread of infection in a population.
The Basic Reproduction Number ℜ0 is obtained by taking the largest (dominant) eigenvalue or spectral radius of ( ) The basic reproduction number is given by (F −1), that is, the highest eigenvalue of Therefore, the basic reproduction number is Theorem 3.2: The disease-free steady state is locally asymptotically stable if 0 1 ℜ < and unstable otherwise.

Basic Idea of Homotopy Analysis Method
In this section, we apply the Homotopy Analysis Method as proposed by Liao (2004), to solve problem (2.1). Consider the following equation: (13) where N is a nonlinear operator, t denotes the independent variable. ( ) is an unknown function. Let 0(t) denote an initial guess of the exact solution ( ), h≠0 an auxiliary parameter, H(t) ≠ 0 an auxiliary function, and L an auxiliary linear operator with the property L [ ( )]= 0 when ( ) = 0. Then using [0,1] as an embedding parameter, we construct such a homotopy.
We have the freedom to choose the initial guess 0(t), the auxiliary linear operator L, the non-zero auxiliary parameter h, and the auxiliary function H(t).
Expanding ∅( ; ) in Taylor series with respect to , we have Where If the auxiliary linear operator, the initial guess, the auxiliary parameter h and the auxiliary function are properly chosen so that a) The solution ∅( ; ) of the zero-order deformation equation (11) According to equation (16), the governing equation can be deduced from the zero-order deformation equation (14).
Define the vector Differentiating the zero-order deformation equation (15) It should be emphasized that 0(t) for m≥1 is governed by linear equation (20) with the linear boundary conditions that comes from the original problem, which can be solved easily symbolic computation software such as Maple, MATLAB or MATHEMATICA.

Solution of the SEIR Model byHomotopy Analysis Method (HAM)
In order to explicitly construct approximate nonperturbative solutions of the system described by equations (1), we employ the concept of homotopy analysis method. The advantage of this method is that it provides a direct scheme for solving the problem. To apply the homotopy analysis method, we choose; S0(t) = Ns, E0(t) = NE, I0(t) = NI and R0(t) = NR.
as initial approximation of S(t), E(t), I(t) and R(t). Let [0,1] be the so-called embedding parameter. The HAM is based on the kind of continuous mappings.
such that as the embedding parameter increases from 0 to 1, ∅1 (t; ) varies from the initial approximation to the exact solution. To ensure this, choose such auxiliary linear operators as By the h-curves, it is reasonable to use hi = −1. Using Hi(t) =1 Hence the solution of the " ℎorder deformation equation of the above equations for m≥1 and using h = −1 and H(t) =1 is given by;

Numerical Solutions
For numerical results, we consider the following values for the parameters involved in the model as presented in Table 2.
Year For high accuracy of results, we use Maple 15 computation software to find the approximate series solution for S(t), E(t), I(t) and R(t), up to the fifth iteration.
The 1st, 2nd, 3rd, 4th and 5th iterations for S(t), E(t), I(t) and R(t) are calculated and the fifth iteration graphed and presented as below; First Iteration

Summary
This work proposed an epidemic model to defend the propagation of an infectious disease. Firstly, we formulated an SEIRS mathematical model and showed that the population classes are non-negative and then we determined the basic reproduction number of the disease by using the next generation matrix method.
Next, with the help of the Basic Reproduction Number, we proved the condition for stability of disease -freeequilibrium.
When the Reproduction Number is less than one, it will cause the disease to die out whereas, if the Reproduction Number is greater than one, this classifies the disease as endemic.
Considering our model and the values stated above in Table 1 the Disease Free Equilibrium is locally asymptotically Stable, which implies that the disease dies out eventually because the Reproduction Number is less than one, hence this was verified by the Simulation results in Figure 3.
We went further to adjust the value of the effective contact rate and got a Reproduction Number greater than one, this implies that Disease Free Equilibrium is locally asymptotically Unstable hence this classifies the disease as Endemic. Also this was verified by the Simulation results in Figure 4.
Finally, it is well known that non-linear Ordinary Differential Equations (ODEs) are difficult to solve than linear ODEs, especially by means of analytic methods.
From our numerical solution, we demonstrated the ability of HAM to converge very fast, we saw that the HAM converges in just few iterations. For the accuracy of the method, it has been established by Liao (1992). Hence, we conclude that HAM is a very efficient and accurate method in solving nonlinear mathematical Model.

Conclusion
After the analysis and simulation study we conclude as follows: The Model has a Disease Free Equilibrium as  Table 1 above we got our Reproduction Number 0 = 0.0081 < 1, hence, on average an infected individual produces less than one new infected individual over the cause of his infectious period and the infection cannot grow , hence the Disease Free Equilibrium (D.F.E.), is locally asymptotically Stable and the disease cannot invade the population. This is also shown in Figure 3 above. In the graph, the Susceptible and Recovered goes up hence there is an increase in the Susceptible and Recovered individuals while the graph of Figure 3 the Exposed and infectious goes down hence there is a decrease in the Exposed an Infectious individuals. Figure 4 above is gotten by changing the effective contact rate I from 0.05 to 7 while all other parameters are left unchanged, hence we got our Reproduction Number 0 = 1.13636 >1 , hence, on average an infected individual produces more than one new infected individual over the cause of his infectious period and the infection can grow, hence, the Disease Free Equilibrium (D.F.E.), is locally asymptotically unstable and the disease can invade the population hence the disease is endemic.
We got our Reproduction Number 0 = 0.0081 < 1, hence this means that on average an infected individual produces less than one new infected individual over the cause of his infectious period and the infection cannot grow , hence the Disease Free Equilibrium (D.F.E.), is locally asymptotically Stable and the disease cannot invade the population. This is also shown in the above, in the graph the Susceptible and Recovered goes up hence there is an increase in the Susceptible and Recovered individuals while thegraph of the Exposed and infectious goes down hence there is a decrease in the Exposed an Infectious individuals.

Recommendation
Globally, diseases have been the major problem to human health, which contributes to the high mortality rate in the world today. Also it thrives in the context of poverty, the costs of seeking accurate diagnosis and treatment can be considerable for low-income household. Infected individuals face substantial cost before diagnosis in that they often consult several public and private providers before and in the process of being diagnosed. Hence in view to these, we therefore recommend that; The Government, non-governmental organizations and stakeholders should help in creating awareness and preventive measures because prevention is better than cure.
The Government should provide more Hospitals in Nigeria to reduce the high population demand in the exiting hospitals.
Infectious Individuals should be properly treated to avoid continuous spread of the disease.
Individuals should shun taking drugs without doctor's prescription, as this will haveanother side effects.
Finally Homotopy Analysis Method is recommended for further work on non-linear Ordinary differential Equations 50.