Dynamics of Toxoplasmosis Disease in Cats population with vaccination

We extend the deterministic model for the dynamics of toxoplasmosis proposed by Arenas et al. in 2010, by separating vaccinated and recovered classes. The model exhibits two equilibrium points, the disease-free and endemic steady states. These points are both locally and globally stable asymptotically when the threshold parameter Rv is less than and greater than unity, respectively. The sensitivity analysis of the model parameters reveals that the vaccination parameter π is more sensitive to changes than any other parameter. Indeed, as expected the numerical simulations reveal that the higher the vaccination rate of susceptible cats the smaller the value of the threshold Rv (i.e., increase in π results in the decrease in Rv), leading to the eradication of toxoplasmosis in cats population. DOI:10.46481/jnsps.2021.141


Introduction
The causative agent of toxoplasmosis is toxoplasma gondii (T. gondii), which is a protozoan of the order Eucoccidioroda and the family Sarcocystidae. T. gondii is characterized as an intracellular parasite that lives in the host cell by regulating vital processes to acquire nutrients, guaranteeing its survival and thus evading the host immune system [12]. It is an infectious pathogen which enters its host through ingestion of either the oocysts, the trachyzoite, or tissue cysts (bradyzoites) from contaminated water, soil, or infected meat. The three major types of disease reservoirs of T. godii are: Domestic and wild cats, soil and water contaminated with feces (Non-living reservoirs) and bradyzoites in tissue cysts of intermediate hosts (animal reservoirs) [25]. The intermediate hosts of T. gondii include sheep, goats, rodents, cattle, swine, chicken and birds. The modes of T. gondii transmissions include ingestion of oocysts that excrete in cats feces for which tissue cysts develop through exposure to cat litter or soil and water. This is the most common and well known modes of transmission to other animals. T. gondii can also be transmitted via consumption of T. gondii tissue cysts in raw or undercooked meats, unpasteurized milk and consumption of oocysts in foods infected by contaminated fomites. This is called the food born transmission. Another route of transmission is congenital toxoplasmos (transplacental transmission) in which a mother infect her offspring in uterus through the blood during pregnancy. A pregnant women with acute infections can transmit to the fetus and cause severe illness such as mental retardation, blindness, and epilepsy [12,25]. It was re-ported in a recent study in dogs that sexual transmission of T. gondii in canine species is possible [25]. Treatment of toxoplasmosis include sulfonamides against murine toxoplasmosis, combined therapy with sulfonamides and pyrimethamine as the standard therapy for toxoplasmosis in humans, spiramycin during pregnancy to reduce transmission of the parasite from mother to fetus, clindamycin especially for patients allergic to sulphonamides [12].
A number of mathematical epidemic models for the transmission dynamics of toxoplasmosis parasites with or without vaccination were proposed in the literature. In 2008, Aranda et al. [6], proposed a simple mathematical model to study the dynamics of Toxoplasmosis disease in the population of Colombia. Numerical simulations of the model using available data reveals the effect of some (hygienic actions, educations programs, more testing and treatments) strategies for the control of toxoplasmosis parasite. This model seems to be the first mathematical model for the transmission dynamics of Toxoplasmosis disease in the human population. In the following year, an epidemic model of toxoplasmosis disease in human and cat populations was presented in [9]. The analysis of the model indicated that toxoplasmosis disease persistence or extinction is determined by the threshold parameter, R 0 . Moreover, they showed by numerical simulations the importance of cats vertical transmission to the dynamics of the infection. In fact, the higher the vertical transmission the higher the proportions of infected cats and infected humans at the endemic state. Thus, the vertical transmission could be an important mechanism that favours the maintenance of the virus areas with low human densities. Another epidemiological model of toxoplasmosis in a cat population with a continuous vaccination schedule was proposed in [2]. Indirectly, the model considers the infection of prey through the oocyst shedding by cats. They proved that the global dynamics and disease outcome are solely determined by the basic reproduction number, R 0 . Furthermore, numerical simulations of the model suggest that the continuous vaccination program is more effective than the removing of oocysts in the environment since the increase of the latter in the environment to feasible values is not enough to reduce the threshold parameter, R 0 to a value less than unity.
The spatial spread of toxoplasmosis parasite was also considered later on. In 2012, a spatial model for the spread of toxoplasma gondii through a heterogeneous predatorprey system was proposed in [16]. They considered some relevant toy models due to the complexity of the proposed model. Then they proved the existence and local stability of a persistent steady state for the underlying predatorprey model systems. A spatial mathematical model for the reproduction dynamics of the Toxoplasma gondii parasite in the definitive host Felis catus (cat) was presented in [10]. Both asexual and sexual reproduction processes of the T. gondii parasite were incorporated in the model. Some numerical results showed that variations in the transition and loss rates do not produce significant changes in the reproduction, propagation and creation of new populations. On the other hand, with either low or high consumption of oocysts from the environment by the cat, the involved populations are always reproduced. Then they spread by all over epithelial cells and subsequently are expelled to the environment through the cat feces. Ferreira et al. [11], studied the dynamical behaviours of both deterministic and spatial models of toxoplasmosis disease in cat and human populations. They showed that the deterministic model exhibits a trans-critical bifurcation and the system has no periodic orbits inside the positive octant. The global asymptotic stability of the endemic equilibrium point of the model was proven in the first octant. Moreover, they carried out the local stability analysis for the spatially homogeneous equilibrium points of the reaction diffusion model. In the restricted region (the first octant), the global stability of both the disease-free and endemic states of such a model were established.
In the last decade, a hybrid model for the spread of toxoplasmosis between two cities of Colombia in cat and human populations was proposed by Peña et al. in [4]. The model was developed using System Dynamics (SD) and Geographic Information Systems (GIS). The analysis of the model indicates that the disease persistence continues in the first city, even after transported to another community. This demonstrated that these diseases should be treated through cooperative mechanisms between communities. A deterministic mathematical model of interaction between Toxoplasma Gondii transmission dynamics and host immune responses was developed and presented in [1]. The local asymptotic stability analysis of the model equilibria and numerical simulations were carried out to understand the transmission dynamics and the impact of different functional responses. In fact, the Holling type II functional response enhances the effect or cells of hosts immune response.
In a recent development, an epidemic mathematical model of toxoplasmosis disease was presented in [15] to study the dynamic behaviour of the parasite for controlling the epidemic. The existence and stability analysis of each equilibrium point of the model was carried out. They constructed a function of two important parameters of the model namely the controlled rate and the rate of infected births which influences the disease spread. The bifurcation analysis of each region was carried out from each function of such parameters. Three regional conditions were obtained through this analysis, showing the dynamics of the toxoplasmosis epidemic of these two important parameters with each interpretation of the bifurcation region. The widespread of co-infection by Plasmodium species and Toxoplasma gondii in humans mostly in the tropical region motivated Ogunmiloro to study the co-infection dynamics of malaria-toxoplasmosis in the mainly human and feline susceptible host population in [18]. The analysis showed that malaria-toxoplasmosis infection free equilibrium is locally and globally stable if the basic reproduction number is less than unity. While the co-infection occurs when the basic reproduction number is greater than unity. Furthermore, sensitivity analysis of the model parameters indicated the need for proper optimum medical strategies to reduce and contain the spread of malaria-toxoplasmosis co-infections.
In 2010, a deterministic mathematical model for the transmission dynamics of toxoplasmosis disease in cats population under vaccination was proposed in [2]. Implicitly, their model considers the infection of prey through the oocyst shedding by 18 cats in the environment. The analysis of the model shows that the global dynamics and disease outcome are completely determined by the basic reproduction number, R 0 . Moreover, the numerical results reveal that the effectiveness of the continuous vaccination program is more than the removing of oocysts in the environment. This is for the fact that an increase of the latter in the environment to feasible values is not enough to control disease persistence by reducing the threshold parameter, R 0 to a value less than one. It is instructive to note that in their model the natural immunity and vaccine protection were assumed to be the same and provide permanent immunity. Therefore, they merged the vaccinated cats and non vaccinated recovered cats in the same class. We extend their model here by separating vaccinated cats from non vaccinated recovered cats to have two classes. This is done for the fact that there is increasing clinical evidence which differentiate between these two immunities against toxoplasmosis in cats [5,7,8,21]. As in [2], we show that the disease outcome, persistence or extinction is determined by the vaccinated reproduction number R v . The sensitivity analysis of the model parameters indicates that the vaccination parameter π is more sensitive to changes than any other parameter of the model. Furthermore, numerical simulations show that the higher vaccination rate of susceptible cats reduces the value of the threshold parameter R v the more. This is also similar to the result presented in [2] that vaccination program is more effective than the removal of oocysts in the environment.
We organise the remaining part of the paper as follows: In Section 2, model description is provided, Section 3 presents the stability analysis of the model equilibria while sensitivity analysis of the model parameters is presented in Section 4. Finally, we discuss the results and give a concluding remark in Section 5.

Model Formulation
To extend the model proposed by Arenas et al. [2], we separate the compartment of vaccinated cats from that of recovered (vaccinated) cats. Thus we divide the total population of cats to four compartments of susceptible (S(t)), cats that may become infected after an effective contact with oocysts, infective (I(t) with T. gondii, vaccinated (V(t)), and non vaccinated recovered R(t). Using the following notations and hypotheses as presented in [2], we have the required model (1). The rates of transfer between the four classes are shown in Figure 1. ii. The birth rate (δ) of cats is assumed to be equal to the cats natural death rate µ, iii. The non vaccinated susceptible cats transits to an infective compartment after an effective contact with oocysts at the rate β. Thus, β depends on the probability of an effective contact with oocysts, vii. µ 0 is the death rate of the oocysts which can be modified by the removal of oocysts from the environment, vi. It is assumed that there is vertical transmission in the infected cats population (i.e., transmission from mother to fetus), viii. Vaccinated cats do not shed oocysts and acquire permanent immunity so that they move to recovered compartment R(t). The recovered cats can be re-infected with T.gondii but do not shed oocysts again and so, they move to susceptible class at the rate τ, ix. The susceptible cats S (t) have the same probability of becoming infected (i.e., homogeneous mixing is assumed) x. κ > 0 is the rate of appearance of new oocysts from the infected cats into the environment, xi. The parameter γ is the progression rate from vaccinated class to the recovered compartment.
In the first equation of model (1), the rate of change of the susceptible cats population is increased by per capita birth rate 19 Rate of appearance of oocyst 1 20 Wu and Pei γ Progression rate of vaccinated cats 0.3 assumed δ of susceptible, vaccinated and recovered classes, the waning immunity after recovery at the rate τ. It decreases by vaccination at the rate, π, natural death µ and effective contact at the rate β between the oocysts and the healthy cats. The rate of change of Infected cats is increased by effective contact of the susceptible cats and the oocysts at the rate β, with birth rate δ of the infected cats. While it decreases by contracting the disease at the rate α. The third equation of model 1 is increased by vaccinating susceptible cats at the rate π and decreases by natural death µ and the rate at which they move to recovered class γ.
Similarly, the rate of change of recovered cats is increased at the per capita rates of recovery of infected and vaccinated cats α and γ, respectively. It is reduced by natural death at the rate µ and wane of immunity from the recovered class at the rate τ. Finally, the rate of change of oocyst with time is increased by shedding the virus at the rate, κ and reduced by the natural death of oocysts at the rate µ 0 . Using R(t) = 1 − S − I − V, the system of equations (1) becomes.
Theorem 1. Model (1) is positively invariant and attractive in the biologically feasible region Proof. Using the first and fourth equations of system (2), we have asymptotically [2]. Hence the region attract all solutions in R 4 + . Thus it is sufficient to consider the dynamics of the model (2) in Ω since it is mathematically well-posed.

Disease Free Equilibrium (DFE)
In the absence of infection (i.e., when all the infected compartments of the model (2) are empty), we have the disease free equilibrium point of the model denoted by P 0 . Then For us to establish the linear stability of equilibria, we need a threshold parameter usually called the basic reproduction number. This threshold parameter is defined as the number of secondary infections produced by one infective when introduce into the susceptible population [13]. The next generation operator method can be used to compute such a threshold parameter as presented in [20]. Moreover, when there is a vaccination such a parameter is called a vaccinated reproduction number (see for instance, [22]). The method consists of two different matrices, the matrix of the new infection terms and that of the transition terms denoted by G and M, respectively. Thus Then Hence, vaccinated reproductive number R v is the largest eigenvalue of the matrix GM −1 . Therefore 3.1.1. Local Stability Analysis of DFE Theorem 2. The DFE, P 0 of model (2) is locally asymptotically stable on Ω if R v < 1 and is unstable whenever R v > 1.
Proof. The Jacobian matrix of system (2) evaluated at the DFE is as follows Then, interchanging row-3 and row-4 and using row operations, the matrix J(P 0 ) becomes Divide row 4 by π, then new row 4 is R 4 = (π + τ + µ)R4 + R1 and so Thus the first two eigenvalues of J(P 0 ) are While the remaining two eigenvalues are for the following 2 × 2 sub-matrix of J(P 0 ).

Global Stability of the DFE
Using Lyapunov Principle method in conjunction with LaSalle's Invariant Principles [19], we establish the global asymptotic stability of the DFE as presented in the following theorem.
Theorem 3. If R v 1, then the DFE, P 0 of model (2) is globally asymptotically stable on Ω.

Endemic Equilibrium (EE)
In the presence of disease infection, the infected compartments (I, O) are non empty and so, the disease invades the population. Let P e = (S e , I e , O e , V e ) be an endemic equilibrium point of model (2), then setting the right hand side of model (2) to zero we obtain after algebraic manipulations.
Hence, P e is the unique endemic equilibrium point of system (2) which exist if and only if R v > 1. 21 3.2.1. Local stability analysis of the EEP Theorem 4. The EEP, P e of model (2) is locally asymptotically stable in the interior of Ω if R v > 1.
Proof. Evaluating the Jacobian matrix of system (2) at P e gives with a = (µ+π+τ+βO e ) and b = βO e . Then we use elementary row operations to obtain Thus the eigenvalues of J(P e ) are Therefore, all the four eigenvalues are real and negative and so, the EEP, P e is locally asymptotically stable if R v > 1. .

Global Stability of Endemic Equilibrium
Theorem 5. If R v > 1, then there exist an endemic equilibrium point P e and it is globally asymptotically stable in the interior of Ω.
Proof. : Given that R v > 1, then the existence of EEP is guaranteed (see Section 3. The time derivative of L along the solution curve of the system (2) yieldṡ But at endemic state, we have µ + τ = (µ + τ)(S e + I e ) + τV e + βS e O e + πS e αI e = βO e S e , πS e = (µ + γ)V e , κI e = µ 0 O e .
Then using (4) in (3), we obtain alter algebraic manipulationṡ Since S ≤ S e , I ≤ I e , V ≤ V e and O ≤ O e , then (5) becomeṡ It follows from arithmetic-geometric inequality that and so,L ≤ 0. Then we conclude this theorem using similar argument as in the proof of Theorem 3.

Sensitivity analysis
One can observe from Section 3, that the disease persistence or otherwise in a community is determined by the threshold parameter R v . Moreover, the uniqueness properties of both disease-free (P 0 ) and endemic (P e ) equilibria indicate that there is an exchange of stability between P 0 and P e when R v = 1 (i.e., P 0 undergoes a transcritical bifurcation at R v = 1, see for instance [23,Theorem 4]). However, it is instructive to investigate the model parameter that has greatest effect on this threshold parameter, thereby having same effect on the toxoplasmosis disease persistence. To achieve that, we can apply sensitivity indices to measure the relative change of a state variable with a changing parameter [17,23].
Then using the definition in [17], the normalised forward sensitivity index with respect to each of the parameters of R v is  Table 1. the following.  Table 1.
It can be observed from equation (7) that most of the expressions for the sensitivity indices are complex and so, we evaluate these indices at the parameter values presented in Table 1. The resulting sensitivity indices of R v to the eight model parameters from the most sensitive to the least one are shown in Table 2.
We infer from Figures 2 and 3 that the threshold parameter R v increases with increasing values of β, κ, µ and τ. While it decreases with increasing values of π, µ 0 and α. As expected the vaccination parameter π is more sensitive to changes than any other parameter of the model.

Discussion and conclusion
In this note, we extend the deterministic model for the transmission dynamics of toxoplasmosis in cats population with vaccination by separating the vaccinated and the recovered compartments. The stability analysis of the model equilibria are carried out and we show that the unique disease-free equilibrium as well as the unique endemic state of the model are asymptotically stable globally under certain condition. That is, the DFE, P 0 is stable globally asymptotically when R v < 1, while the EEP, P e is globally asymptotically stable if R v > 1. In order to complete our investigation, we use the sensitivity indices as presented in [17] to determine the most sensitive model parameter. As expected, the vaccination parameter, π is the most sensitive to changes than any other parameter of the model. Moreover, considering the complex nature of the expressions for the sensitivity indices in equation (7), we evaluate these indices at the parameter values presented in Table 1. Then, we order them from the most sensitive to the least one as presented in Table 2.
Furthermore, the numerical simulations in Figures 2 and 3 reveal how the threshold parameter R v varies with varying values of each of the model parameters. In fact, R v increases with increasing values of the parameters β, κ, µ and τ and decreases with increasing values of π, µ 0 and α.
In conclusion, the sensitivity analysis indicates that the best way to contain the spread of toxoplasmosis is either to increase the vaccination rate of susceptible cats or to increase the death rate of the oocysts by environmental sanitation. It can be observed from Figure 2 that vaccination is the best control strategy than the removal of the oocysts. This is similar to the result presented in [2].