J. Nig. Soc. Phys. Sci. 4 (2022) 193–204 Journal of the Nigerian Society of Physical Sciences Mathematical Model of In-host Dynamics of Snakebite Envenoming S. A. Abdullahia,c, A. G. Habibb, N. Hussainic,∗ aDepartment of Mathematics, Modibbo Adama University Yola, Adamawa State, Nigeria bInfectious and Tropical Diseases Unit, Department of Medicine, Bayero Univesrity Kano, Nigeria cDepartment of Mathematical Sciences, Bayero University Kano, P.M.B. 3011, Kano, Nigeria Abstract In this paper, we develop an in-host mathematical model of snakebite envenoming that includes tissue, red blood and platelet cells of humans as specific targets of different kinds of toxins in the snake venom. The model is use to study some harmful effects of cytotoxic and hemotoxic snake venom on their target cells under the influence of snake antivenom. The model has two equilibrium points, namely, trivial and venom free. It has been shown that both the equilibrium points are globally asymptotically stable and numerical simulations illustrate the global asymptotic stability of the venom free equilibrium point. Furthermore, simulations reveal the importance of administering antivenom to avert the possible damage from venom toxins on the target cells. It is also shown through simulation that administering the required dose of antivenom can lead to the elimination of venom toxins within one week. Therefore, we recommend the administration of an adequate dose of antivenom therapy as it helps in deactivating venom toxins faster and consequently enhances the recovery time. DOI:10.46481/jnsps.2022.548 Keywords: Snakebite, in-host model, Venom, Antivenom, Stability analysis Article History : Received: 21 December 2021 Received in revised form: 17 February 2022 Accepted for publication: 25 February 2022 Published: 29 May 2022 c©2022 Journal of the Nigerian Society of Physical Sciences. All rights reserved. Communicated by: T. Latunde 1. Introduction Snakebite envenoming is one of the most devastating neglected tropical diseases. It continues to impose a major public health and socio-economic burden in low and middle income coun- tries around the world [1, 2, 3, 4, 5]. This neglected tropical disease, which is considered a serious public health issue in tropical and subtropical countries worldwide, is caused by ven- omous snake bites [6, 7, 8]. Globally, more than five million people are bitten by snakes every year [3, 4, 9]. It has been estimated that 1.8 - 2.7 million envenomations and 81 000 - ∗Corresponding author tel. no: +2348037594137 Email address: nhussaini.mth@buk.edu.ng (N. Hussaini) 138 000 deaths occur annually, with an additional 400 000 am- putations and other severe health consequences, such as infec- tion, tetanus, scarring, contractures, and psychological sequelae [3, 4, 6, 7, 8, 10, 11]. The majority of snakebites occur in Africa and Southeast Asia. Snakebites are most common among peo- ple living in rural, resource-poor settings. Furthermore, agricul- tural workers, women, and children are the groups most at risk of being bitten by snakes [6, 7, 8, 9, 12, 13]. According to the World Health Organization (WHO) [14], there are more than 3000 species of snakes in the world, of which 600 are venomous and more than 200 are medically important. Venomous snakes cause substantial morbidity and mortality in humans. When a venomous snake strikes, it may likely inject venom into its vic- 193 Abdullahi et al. / J. Nig. Soc. Phys. Sci. 4 (2022) 193–204 194 tim, as some bites may be dry, which is a bite that will not lead to the release of venom into its target [15]. These venomous snakes administer bites either as a method of hunting their prey, or as a means of protection against their predators [3]. Accord- ing to Young and Zahn in [16] the rate of venom discharge and total volume of venom ejected by venomous snakes is much greater during defensive strikes than during predatory strikes. A defensive strike can have 10 times as much venom volume expelled at 8.5 times the flow rate. Snake species accountable for the majority of envenoming in the human population are categorized within the families of Elapidae (cobra, king cobra, krait, etc.) and Viperidae (Saw-scaled viper, pit viper, Russell’s viper). Further, species belonging to the families Atractaspidi- dae and Colubridae can also induce significant envenomation [17, 18, 19, 20, 21, 22, 23]. Snake venom is the poisonous, nat- urally yellow fluid stored in the adapted salivary glands of ven- omous snakes. Snake venom is mainly categorized into three namely: cytotoxins, neurotoxins and hemotoxins, which tar- get body cells, nervous system and cardiovascular system re- spectively [14, 24]. According to León et al. in [25] when a snake bites and injects venom into its victim, two concurrent processes occur within the organism: (a) the development of toxic effects and (b) stimulation of immune responses in order to neutralize venom proteins. If the capability of snake venom to cause local and systemic mutilation surpasses the ability of the animal to assimilate and respond to the violence, an enven- omation is established. According to some studies, the major- ity of snakebites occur on the hands, arms, or legs [26, 27]. Some signs and symptoms of snakebite may include: redness, swelling, and severe pain at the bite site, which may take up to an hour to appear. Further, vomiting, blurred vision, tin- gling of the limbs, and sweating may likely occur [26, 28]. The severity of snakebite depends on some factors, such as the kind of snake, the part of the body bitten, the quantity of venom injected, and the general health of the person bitten. The severity of snakebite is higher in children than in adults, due to their smaller size [9, 13, 29]. Snakebites can be avoided by wearing protective footwear, avoiding snake-infested areas, and not handling snakes [7, 13, 28, 30]. The use of snake an- tivenom remains the most effective method of averting death from bites. Nevertheless, antivenoms often have side effects [9, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48]. Some researches on snake venom focused on the use of venom toxins to control infections such as protozoal infections, cancers, and so on (see for instance [49, 50, 51, 52, 53] and references therein). Numerous population-based mathematical models have been developed to study the transmission dynam- ics and control of neglected tropical diseases and other diseases (see, for instance, [54, 55, 56, 57, 58, 59, 60] and references therein). Also, several in-host mathematical models for dis- eases such as HIV, HBV, HCV, Chagas disease, etc. have been successfully developed and were used to study the transmis- sion dynamics and possible control of the pathogens inside the victim’s body (for instance, see [61, 62, 63, 64, 65, 66] and ref- erences therein). However, the literature on the use of mathe- matical modeling techniques to gain an in-depth understanding of the dynamics of snake venom inside an envenomed person is scarce. Although not a mathematical model, Kini and Evans in [67] proposed a hypothetical model and explained how venom phospholipase A2 enzymes exhibit specific pharmacological ef- fects such as presynaptic neurotoxicity, cardiotoxicity, myotox- icity, anticoagulant and platelet effects. The research revealed how toxins and enzymes contained within the venom are drawn to specific tissues or cells due to their affinity for specific pro- teins. In 1994, Stagg [68] developed a mathematical model and studied the effect of cobra venom and its anti-venom inside the envenomed body (in-vivo). The work primarily described the neurotoxic type of snake venom and the cell that it targets. Furthermore, Tanos et al. [69] developed a semi-mechanistic model of the clotting cascade to describe the effect of procoag- ulant toxin from taipan venom and the effect of anti-venom in controlling venom-induced consumptive coagulapathy. Their work focuses on how the hemotoxic venom affects or disrupts blood clotting factors. We are motivated by the fact that in June 2017, WHO added snakebite envenoming to its priority list of neglected tropical diseases and has set a target to reduce its mor- tality and disability by 50% before the year 2030 [70]. Also, to the best of our knowledge, not much work has been done in terms of mathematical modeling to gain more insight into the dynamics and control of toxic effects of snake venom inside the host body. Thus, the aim of this study is to develop a novel in- host mathematical model that will describe the dynamics and control of cytotoxic snake venom that was not considered by the aforementioned studies. Moreover, the study will also in- vestigate the impact and control of hemotoxic snake venom on red blood cells and platelets. Additionally, it is important to note that providing effective control of snakebite envenoming in a population may depend on a comprehensive understanding of the effects and dynamics of various kinds of toxins inside snake venom on their target cells (in-vivo). Therefore, this pa- per is committed to study the effects of interactions between snake venom (cytotoxic and hemotoxic) and their associated target cells with and without antivenom therapy. The paper is organized as follows: the model is formulated in section 2, analyzed in section 3. Results and discussion are reported in section 4 and conclusion is provided in section 5. 2. Model Formulation A mathematical model that consists of nine in-host compart- ments is proposed in order to investigate the effects of snake venom and its target cells inside the human victim at a given time. The in-host compartments are: quantity of healthy tissue cells denoted by T (t), quantity of healthy red blood cells rep- resented by R(t), quantity of healthy platelets denoted by P(t), quantity of snake venom with cytotoxic and hemotoxic com- ponents expressed as V (t), quantity of antivenom to be admin- istered denoted by A(t). We assumed that the quantity of an- tivenom to be administered is sufficient to neutralize the circu- lating venom within the victim’s bloodstream. We also included damaged tissue cells (necrotic cells) caused by the cytotoxic component of venom toxins (necrotoxins), denoted by N(t), damaged red blood cells (hemolyzed cells) caused by the hemo- toxic component of venom toxins, denoted by H(t), damaged 194 Abdullahi et al. / J. Nig. Soc. Phys. Sci. 4 (2022) 193–204 195 platelets caused by the hemotoxic component of venom toxins (platelets toxins), denoted by D(t), and venom-antivenom com- plexes expressed as C(t). The quantity of healthy tissue cells is replenished at a logistic growth rate given by γT T ( 1 − TKT ) with T ≤ KT . It is decreased as it undergoes necrosis by direct action of the venom’s cytotoxic components (necrotoxin) that begin at the bite site at a constant destruction rate that follows the law of mass action given by (β1T V ). The damaged cells caused by venom-induced necrosis are known as necrotic cells. The healthy tissue cells further decreased due to the elimination of the dead cells at a rate µT . Thus, the equation governing the dynamics of healthy tissue cells is given by dT dt = γT T ( 1 − T KT ) −β1T V −µT T. (1) The quantity of healthy red blood cells is replenished at a lo- gistic growth rate denoted by γRR ( 1 − RKR ) (R ≤ KR). After 120 days it is reduced by natural elimination from the host body at a rate µR. It further decreased because of the damage it undergoes as a result of the deleterious effects of the hemotoxic component of the venom as it circulates in the victim’s blood stream. The action of the hemotoxic component of the venom also follows the law of mass action given by (β2RV ). The damaged red blood cells due to venom-induced hemolysis are termed hemolyzed cells. Thus, we have dR dt = γRR ( 1 − R KR ) −β2RV −µRR. (2) The quantity of healthy platelets is replenished at a logistic growth rate of γP P ( 1 − PKP ) (P ≤ KP). It is reduced by natu- ral elimination of the dead platelet cells at a rate of µP. It de- creases further because the hemotoxic component of the venom (platelet toxins) activates and damages platelets at a rate of β3. This effect may cause serious internal bleeding. Thus, the fol- lowing equation is obtained dP dt = γP P ( 1 − P KP ) −β3 PV −µP P. (3) It is assumed that an initial quantity of venom injected by an offending snake is in circulation in the host body. Therefore, the quantity of venom decreases because of its reaction with the various target cells, i.e. tissue, red blood, and platelet cells at the rates β1,β2 and β3 respectively. It further diminishes due to its neutralization by antivenom that is describe by the law of mass action (β4 AV ). This reaction produces venom-antivenom complexes that are not harmful to the body. The amount of venom is also eliminated at a rate µV . Therefor, the following equation is formed dV dt = −β1T V −β2RV −β3 PV −β4 AV −µV V. (4) When the required dose of snake antivenom is administered, its action is to neutralize the circulating venom toxins inside the victim’s body. The quantity of antivenom inside the host body is decreased because of the reaction between venom and antivenom that follows the law of mass action given by (β4 AV ). The unused antivenom is then removed at a rate of µA. Thus, we have dA dt = −β4 AV −µA A. (5) The quantity of necrotic cells is formed because of the reaction between the cytotoxic component of the venom (necrotoxin) and the tissue cells at the rate β1 and is eliminated at a rate µN . Therefore, we obtain dN dt = β1T V −µN N. (6) The quantity of hemolyzed cells is formed due to the reaction between the hemotoxic component of the venom and the red blood cells at the rate β2 and is eliminated at a rate µH . Thus, we have dH dt = β2RV −µH H. (7) The formation of the quantity of damaged platelets occurs as a result of reaction between the platelet toxins and the healthy platelet cells at the rate β3. It decays at a rate of µD. Therefore, we have dD dt = β3 PV −µD D. (8) The quantity of venom-antivenom complexes is produced be- cause of the neutralization of venom by antivenom and is elim- inated from the body at a rate µC . So that dC dt = β4 AV −µCC. (9) Based on the aforementioned description of the model, the in- host dynamics of interaction between cytotoxic and hemotoxic venom, their target cells, and antivenom is described by the following system of first order nonlinear ordinary differential equations. The schematic diagram is depicted in Figure 1, and the corresponding variables and parameters are tabulated in Ta- bles 1 and 2, respectively: dT dt = γT T ( 1 − T KT ) −β1T V −µT T, dR dt = γRR ( 1 − R KR ) −β2RV −µRR, dP dt = γP P ( 1 − P KP ) −β3 PV −µP P, dV dt = −β1T V −β2RV −β3 PV −β4 AV −µV V, dA dt = −β4 AV −µA A, dN dt = β1T V −µN N, dH dt = β2RV −µH H, dD dt = β3 PV −µD D, dC dt = β4 AV −µCC, (10) 195 Abdullahi et al. / J. Nig. Soc. Phys. Sci. 4 (2022) 193–204 196 Figure 1. Schematic diagram of the model. Note that, in the schematic diagram of the model, the dotted lines indicate the reactions between venom toxins and their target cells (which produce damaged cells) and the neutralization of venom tox- ins by antivenom (which leads to the formation of venom-antivenom complexes). The solid arrows represent the self-replenishment of the healthy cells, the formation of damaged cells and venom-antivenom complexes. It also indicates the natural elimination of damaged or dead cells from the host body. Table 1. Description of State Variables of the Model. Variables Description T(t) Quantity of healthy tissue cells at time t R(t) Quantity of healthy red blood cells at time t P(t) Quantity of healthy platelets at time t V(t) Quantity of venom in the victim blood stream at time t A(t) Quantity of antivenom in the victim blood stream at time t N(t) Quantity of damaged tissue (necrotic cells) at time t H(t) Quantity of damaged red blood cells (hemolyzed cells) at time t D(t) Quantity of damaged platelets at time t C(t) Quantity of venom-antivenom complexes at time t with the initial conditions T (0) > 0, R(0) > 0, P(0) > 0, V (0) = V0 ≥ 0, A(0) = A0 ≥ 0, N(0) ≥ 0, H(0) ≥ 0, D(0) ≥ 0, C(0) ≥ 0. (11) 3. Model Analysis 3.1. Basic properties of the model This section presents some essential properties of the in-host model. It is vital to show that the state variables of the model Table 2. Description of Parameter of the Model. Parameter Description γi(i = T, R, P) Replenishment rates of tissue cells, red blood cells and platelets respectively Ki(i = T, R, P) Maximum quantities of tissue cells, red blood cells and platelets respectively µi(i = T, R, P) Elimination rates of tissue, red blood, platelet cells respectively µi(i = N, H, D) Elimination rates of damaged tissue, red blood, and platelet cells respectively µi(i = V, A, C) Elimination rates of venom, antivenom and venom-antivenom complexes respectively β1 Rates at which cytotoxic component of the venom (necrotoxins) destroys tissue cells β2 Rate at which hemotoxic component of the venom (hemotoxins) destroys red blood cells β3 Rate at which hemotoxic component of the venom (platelet toxins) damages platelets β4 Naturalization rate of venom by antivenom are non-negative for all t > 0. If this result is established, then it is sufficient to say that the model (10) is mathemati- cally and biologically reasonable. It is crucial to assume that the model parameters given in Table 2 are positive. For math- ematical convenience let µ1 = min{µT ,µN}, µ2 = min{µR,µH} and µ3 = min{µP,µD}. Further, let V0 denotes the amount of venom that the offending snake injected into its victim at the initial time. Thus, we assume that quantity of venom inside the host at any given time is presented as V (t) ≤ V0. (12) Similarly, let A0 be the initial dose of antivenom administered into the victim’s bloodstream, which we assumed was sufficient to neutralize the circulating venom for simplicity. Thus, it is assumed that the quantity of antivenom inside the host at any given time is given A(t) ≤ A0. (13) Theorem 3.1. Let the initial data of the model (10) be non- negative, then the model solutions given by (T, R, P, V, A, N, H, D, C) are all non-negative for every t > 0. Proof. Let t1 = sup{t > 0 : T > 0, R > 0, P > 0, V > 0, A > 0, N > 0, H > 0, D > 0, C > 0, }. Thus, t1 > 0. Taking the first equation of the model (10) we have, dT dt = γT T ( 1 − T KT ) −β1T V −µT T. (14) Since, T ≤ KT we have dT dt ≥−(β1T V + µT T ) . (15) 196 Abdullahi et al. / J. Nig. Soc. Phys. Sci. 4 (2022) 193–204 197 (a) (b) (c) Figure 2. Graph showing the lethal effects of Russell’s Viper venom with and without antivenom on (a) tissue cells (b) red blood cells (c) Platelets. In Figures (a), (b) and (c), the red dotted lines indicate the quantities of damaged tissue cells, red blood cells and platelets without the influence of antivenom. While the blue solid lines represent quantities of damaged tissue cells, red blood cells and platelets under the influence of antivenom. Now integrating both sides from t=0 to t = t1 we obtain,∫ t1 0 ( 1 T ) dT ≥− ∫ t1 0 (µ) dt + ∫ t1 0 (β1V (τ)) dτ  T (t) ≥ T (0) exp − µt1 + ∫ t1 0 β1V (τ)dτ   > 0. Thus, T ( t1 ) > 0 for all t 1 > 0, hence, T (t) > 0 for t > 0. Following similar technique, it can be established that the re- maining eight variables R, P, V, A, N, H, D, C are all positive for t > 0. Consequently , the solutions of the system (10) remain positive for all t > 0. Lemma 3.1. The closed set Ω = {(T, R, P, V, A, N, H, D, C) ∈ R9+ : P1 ≤ γT KT µT , P2 ≤ γR KR µR , P3 ≤ γP KP µP , C ≤ β4 A0 V0 µC }, is positively invariant. Proof. Let P1 represents the total quantity of healthy and dam- aged tissue cells (necrotic). Thus, P1 = T + N and dP1 dt = dT dt + dN dt . Now, adding the first and sixth equations of system (10) we get dP1 dt = γT T ( 1 − T KT ) −µ1 P1. (16) It follows that dP1 dt ≤ γT KT −µ1 P1. Applying comparison theorem [73] it can be shown that P1(t) ≤ P1(0)e −µ1 t + γT KT µ1 [ 1 − e−µ1 t ] . 197 Abdullahi et al. / J. Nig. Soc. Phys. Sci. 4 (2022) 193–204 198 (a) (b) (c) Figure 3. Simulations showing the effects of varying the initial amount of Russell’s Viper venom injected on (a) tissue cells (b) red blood cells (c) Platelets. The red dashed, black solid, green dashed, and blue dashed lines denote the quantities of damaged tissue, red blood and platelet cells when 0mg, 50mg, 100mg and 150mg of venom are injected, respectively. Thus, lim sup t→∞ P1(t) ≤ γT KT µ1 . Similarly, let P2 denotes the total quantity of healthy and dam- aged red blood cells (hemolyzed). Thus, P2 = R + H and dP2 dt = dR dt + dH dt . Now, adding the second and seventh equations of system (10) we obtain dP2 dt = γRR ( 1 − R KR ) −µ2 P2. (17) It follows that dP2 dt ≤ γR KR −µ2 P2. Applying comparison theorem [73] it can be shown that P2(t) ≤ P2(0)e −µ2 t + γR KR µ2 [ 1 − e−µ2 t ] . Thus, lim sup t→∞ P2(t) ≤ γR KR µ2 . In the same manner let P3 denotes the total quantity of healthy and damaged platelets. Thus, P3 = P + D and dP3dt = dP dt + dD dt . Now, adding the third and eighth equations of system (10) we get dP3 dt = γP P ( 1 − P KP ) −µ3 P3. (18) It follows that dP3 dt ≤ γP KP −µ3 P3. Applying comparison theorem [73] it can be shown that P3(t) ≤ P3(0)e −µ3 t + γP KP µ3 [ 1 − e−µ3 t ] . Thus, lim sup t→∞ P3(t) ≤ γP KP µ3 . Finally, taking the ninth equation of 198 Abdullahi et al. / J. Nig. Soc. Phys. Sci. 4 (2022) 193–204 199 (a) (b) (c) (d) Figure 4. Simulations showing the effects of varying the initial dose of snake antivenom in averting damaged (a) tissue cells (b) red blood cells (c) Platelets and (d) deactivating venom inside the host body. In Figures 4a, 4b and 4c, the red dashed, black solid, green dashed, and blue dashed lines represent the quantities of necrotic, hemolyzed and damaged platelet cells when 0ml, 100ml, 150ml and 250ml of antivenom are administered, respectively. While in Figure 4d, the black dashed, blue dashed and red solid lines represent the quantities of venom inside the host body when 50ml, 100ml and 250ml of antivenom are administered, respectively. system (10) and using V (t) ≤ V0, A(t) ≤ A0 we have dC dt ≤ β4 A0V0 −µCC. (19) Consequently, using the comparison theorem in conjunction with the integrating factor technique of solving first order linear differential equations, we have C(t) ≤ C(0)e−µC t + β4 A0V0 µC [ 1 − e−µC t ] . Thus, lim sup t→∞ C(t) ≤ β4 A0 V0 µC . Hence, Ω is positively invariant, meaning that, all the solutions starting in Ω stay in Ω for t > 0. Thus, it is established that Ω is an attracting set. Since, the region Ω is positively invariant, it suffice to inves- tigate the dynamics of the model equations given by (10) in Ω, where the usual existence, uniqueness, continuation of solu- tions hold. 4. Results and Discussion 4.1. Equilibrium Points of the Model To obtain the equilibrium points of the model equations given by system (10) the right hand side of each equation in system (10), is set to zero and solves for the associated state variables. The model equations have two non-negative equilibrium points, viz; the trivial and the venom free equilibrium points, denoted 199 Abdullahi et al. / J. Nig. Soc. Phys. Sci. 4 (2022) 193–204 200 by ε0 and ε1 correspondingly. The results obtained are pre- sented in the following subsections. 4.1.1. Trivial Equilibrium Point The trivial equilibrium point is presented as follows: ε0 = ( T 0, R0, P0, V 0, A0, N0, H0, D0, C0 ) = (0, 0, 0, 0, 0, 0, 0, 0, 0) . It is important to note that this equilibrium point is biologically not feasible. However, its qualitative analysis is of interest. 4.1.2. Venom Free Equilibrium Point The venom free equilibrium point of the model is given by: ε1 = (T ∗, R∗, P∗, V∗, A∗, N∗, H∗, D∗, C∗) = ( KT (γT −µT ) γT , KR (γR −µR) γR , KP (γP −µP) γP , 0, 0, 0, 0, 0, 0 ) . (20) This equilibrium point is biologically and mathematical feasi- ble. It describes the healthy state of tissue, red blood cells and platelets. Observe that if γT = µT ,γR = µR and γP = µP then the equilibrium point ε1 reduces to ε0. Therefore, for the equi- librium point ε1 to be positive we must have γT > µT ,γR > µR and γP > µP. 4.1.3. Global Stability of Equilibrium Points of the Model Here we shall use Lyapunov function theory in conjunction with La-Salle’s invariance principle to establish the global asymp- totic stability of the trivial and venom free equilibrium points. 4.1.4. Global Stability of Trivial Equilibrium Point Theorem 4.1. The trivial equilibrium point (ε0) of model (10), is globally asymptotically stable in Ω. Proof. The stability analysis of (ε0) is simplified by defining a perturbation W = X −ε0, X = (T, R, P, V, A, N, H, D, C). (21) Therefore, model (10) is expressed as: dW dt = B(W)W, (22) where W = W Ti (i = 1, ..., 9) and B(W) is a matrix of coefficients of the system (10) with the variables Wi(i = 1, ..., 9), this is pre- sented as follows: where, Z11 = (γT−β1 W4−µT )KT−2 γT W1 KT , Z22 = (γR−β2 W4−µR )KR−2 γR R KR , Z33 = (γP−β3 W4−µP )KP−2 γP P KP , Z44 = −β1W1 − β2W2 − β3W3 − β4W5 −µV. Observe that εW = (0, 0, 0, 0, 0, 0, 0, 0, 0) , is the only equilib- rium point of equation (22). Furthermore, we consider the fol- lowing Lypunov function as applied in the following studies [74, 75, 76] and references therein. V (W) = 〈Y, W〉, Y = (1, 1, 1, 1, 1, 1, 1, 1, 1) > 0. (24) V′(W) = 〈Y, B(W)W〉, = Z11W1 + Z22W2 + Z33W3 − (W1β1 + W2β2 + W3β3 + W5β4 + µV ) W4 − (W4β4 + µA) W5 −µN W6 −µH W7 −µDW8 −µC W9, = ( γT ( 1 − W1 KT ) − γT W1 KT −β1W4 −µT ) W1 + ( γR ( 1 − W2 KR ) − γRW2 KR −β2W4 −µR ) W2 + ( γP ( 1 − W3 KP ) − γPW3 KP −β3W4 −µP ) W3 − (β1W1 + β2W2 + β3W3 + β4W5 + µV ) W4 − (β4W4 + µA) W5 −µN W6 −µH W7 −µDW8 −µDW9. (25) Consider the following conditions: W1 ≤ KT , W2 ≤ KR, W3 ≤ KP. (26) Therefore, applying the conditions given in (26) in equation (25) we obtain V′(W) ≤− ( γT W1 KT + β1W4 + µT ) W1− ( γRW2 KR + β2W4 + µR ) W2 − ( γPW3 KP + β3W4 + µP ) W3 − (β1W1 + β2W2 + β3W3 + β4W5 + µV ) W4 − (β4W4 + µA) W5 −µN W6 −µH W7 −µDW8 −µC W9. (27) Thus, V′(W) ≤ 0 with V′(W) = 0, if and only if Wi = 0, (i = 1, 2, ..., 9). Hence, V (W) is Lypunov function in Ω. Fur- ther, since, Ω is an invariant and attracting set as established in Lemma (3.1) it follows from La-Salle’s invariance principle in [77] that the maximum invariant set contained in {V/V′(W) = 0} is εW . Consequently, the transformed equilibrium point εW, is globally asymptotically stable in Ω. Hence, the trivial equilib- rium, ε0, is also globally asymptotically stable in Ω. 4.1.5. Global Stability of Venom Free Equilibrium Point Theorem 4.2. The venom free equilibrium point (ε1) of model (10), is globally asymptotically stable in Ω. Proof. Consider the following linear Lyapunov function given by: W (t) = a1V (t) + a2 N (t) + a3 H (t) + a4 D (t) + a5C (t) ,(28) where ai, i = 1, ...5 are positive constants to be chosen later. The time derivative of the Lyapunov function along the solutions of model (10) is given as; Ẇ (t) = a1V̇ (t) + a2 Ṅ (t) + a3 Ḣ (t) + a4 Ḋ (t) + a5Ċ (t) , = a1 [ −β1T V −β2RV −β3 AV −µV V ] + a2 [ β1T V −µN N ] + a3 [ β2RV −µH H ] + a4 [ β3 PV −µD D ] + a5 [ β4 AV −µCC ] . (29) 200 Abdullahi et al. / J. Nig. Soc. Phys. Sci. 4 (2022) 193–204 201 B (W) =  Z11 0 0 −β1W1 0 0 0 0 0 0 Z22 0 −β2W2 0 0 0 0 0 0 0 Z33 −β3W3 0 0 0 0 0 −β1W4 −β2W4 −β3W4 Z44 −β4W4 0 0 0 0 0 0 0 −β4W5 −β4W4 −µA 0 0 0 0 β1W4 0 0 β1W1 0 −µN 0 0 0 0 β2W4 0 β2W2 0 0 −µH 0 0 0 0 β3W4 β3W3 0 0 0 −µD 0 0 0 0 β4W5 β4W4 0 0 0 −µC  , (23) Ẇ = (a1 − a2) β1T V + (a1 − a3) β2RV + (a1 − a4) β3 PV + (a1 − a5) β4 AV − a1µV V − a2µN N − a3µH H − a4µD D − a5µCC. (30) Choosing ai = 1 µV µNµHµDµC for (i = 1, ...5) and applying it in equation (30) we obtain Ẇ = − [( 1 µNµHµDµC ) V + ( 1 µVµHµDµC ) N + ( 1 µVµNµDµC ) H ] − [( 1 µVµNµHµC ) D + ( 1 µVµNµHµD ) C ] . (31) Consequently, using equation (31), it follows that Ẇ < 0 un- conditionally, since all the model parameters are positive in Ω. For Ẇ = 0 we must have V = N = H = D = C = 0. Following the claim that Ω is invariant and attracting set, it follows that the largest possible invariant set in Ω is the singleton ε1. According to La-Salle’s invariance principle [77], ε1 is globally attractive. Therefore, the venom free equilibrium, ε1, is globally asymp- totically stable. 4.2. Numerical Simulation In this section, some numerical simulations were performed us- ing the values of the model parameters in Table 3. According to Bawaskar in [79], the deadly dose of Russell’s viper venom is V0 = 150mg and the initial dose of antivenom needed to neu- tralize the venom is A0 = 250ml. The numerical simulation results are shown in Figures 2, 3, 4, 5 and 6. 4.2.1. Effects of venom toxins on target cells with and without antivenom Figure 2a shows the deleterious effects of toxins inside the Rus- sell’s viper venom on the tissue cells. It can be observed that the quantity of damaged tissue cells is on the increase when an- tivenom is not administered (red dot line). On the other hand, a negligible quantity of damaged tissue cells is noticed when an- tivenom is administered (blue line). Also, in Figures 2b and 2c Table 3. Values of parameter used in numerical simulation. Parameter Values Reference γT 0.009 week −1 Estimated γR 1 17.143 week −1 [71, 72] γP 1 1.429 week −1 [78] Ki(i = T, R, P) 80000, 10000, 5000 Estimated µT 1 3 week −1 Estimated µR 1 17.143 week −1 [71, 72] µP 1 1.429 week −1 [78] µi(i = V, A, N) 0.1, 0.013, 0.001 Estimated µi(i = H, D, C) 0.0012, 0.00013, 0.0002 Estimated βi(i = 1, 2) 7.1 × 10−11, 5.1 × 10−06 Estimated βi(i = 3, 4) 3.7 × 10−04, 0.87 Estimated we notice a significant amount of damaged red blood cells and platelets in the absence of antivenom to neutralize the circulat- ing venom. However, an insignificant number of damaged red blood cells and platelets are observed when snake antivenom is administered. 4.2.2. Effects of varying initial quantity of venom on the target cells The dynamical effect of increasing the initial dose of venom on the target cells is depicted in Figure 3. The outcome reveals that the amounts of injured tissue, red blood cells, and platelet cells increase as the initial dose of venom increases. This re- sult suggests that the magnitude of damage to be inflicted by venom toxins on the target cells depends on the initial quantity of venom injected by the offending snake. Therefore, deter- mining the amount of venom injected by the offending snake is crucial in the management of snakebite patients. 4.2.3. Effects of varying initial dose of antivenom on the dam- aged cells and the venom The simulations in Figures 4a, 4b and 4c demonstrate the effect of varying antivenom doses in preventing damaged cells as a re- sult of venom toxins. It is obvious that the application of the re- 201 Abdullahi et al. / J. Nig. Soc. Phys. Sci. 4 (2022) 193–204 202 Figure 5. Simulation showing the effect of using adequate dose of an- tivenom in decreasing the quantity of venom inside the host. quired dose of antivenom is more effective in terms of averting damaged cells. Also, the dynamics of varying the initial dose of antivenom in deactivating venom toxins inside the host body is depicted in Figure 4d. It is evident that when a lower-than- required amount of antivenom is administered, it takes longer time to neutralize the quantity of venom inside the host. On the other hand, when the required amount of antivenom is ap- plied, it takes less time to neutralize the venom toxins. This finding supports the need for administering an adequate dose of antivenom in order to enhance the time of neutralizing the lethal dose of venom within the host. 4.2.4. Effect of adequate dose of antivenom in decreasing the quantity of venom inside the host body According to the simulation result shown in Figure 5, applica- tion of adequate dose of antivenom leads to the elimination of a deadly dose of Russell viper venom within one week. This result is in line with the recovery time for snakebite envenom- ing in Treatment and Research Hospital, Katungo, Gombe state, as reported in [80]. This finding supports the use of the recom- mended dose of antivenom in the treatment of snakebites enven- oming, since it will not only prevent harm to the cells but also speed up recovery time. It is worth mentioning that the epi- demiological implication of the numerical results provided in Figures 2, 3, 4 and 5 is that whenever snakebite envenomation is confirmed, the required doses of antivenom therapy should be administered as soon as possible to avoid the potential dam- age that venom toxins will inflict on target cells. As a result, envenomation-related deaths and disabilities will be avoided. 4.2.5. Numerical illustration of global asymptotic stability of the venom free equilibrium point The results in Figures 6a, 6b and 6c illustrate the convergence of solutions to the venom free equilibrium point which uphold the result of global asymptotic stability of venom free equilibrium point presented in Theorem 4.2. (a) (b) (c) Figure 6. Simulation results showing total quantity of damaged (a) tis- sue cells (b) red blood cells (c) platelets. In all the cases different initial conditions were used to illustrate the result of global asymptotic stability of the venom free equilibrium point. 5. Conclusions In this work, we developed and analyzed an in-host model of snakebite envenoming. The model was used to study some toxic effects of snake venom on tissue, red blood, and platelet cells of humans under the influence of snake antivenom. The basic properties of the model in terms of positvity and bound- edness of solutions were explored. The qualitative study of the model reveals that the model has two non-negative equilibrium points, namely, trivial and venom- free. The equilibrium points were shown to be globally asymptotically stable and numeri- cal simulations were used to illustrate the result of the global asymptotic stability of the venom free equilibrium point as es- tablished in Theorem 4.2. Furthermore, numerical simulations also revealed the importance of antivenom therapy in prevent- ing the deleterious effects of snake venom on target cells. It is also shown through simulation that the application of re- quired dose of antivenom can lead to the elimination of venom toxins within one week. The public health implication of the simulation results in Figures 2, 3, 4 and 5 in terms of manag- ing snakebite envenomation is that determining the quantity of venom injected by the offending snake is crucial. In addition, the required dose of antivenom should be administered as soon as possible so that envenomation-related deaths and disabili- ties can be prevented. Some of the numerical findings in this study provide some useful insights into the management and control of snakebite envenoming patients, such as determining the initial quantity of venom injected by the offending snake and applying the required dose of antivenom to neutralize the venom toxins. This work is limited to the effect of cytotoxic and hemotoxic venoms on the tissue, red blood, and platelet 202 Abdullahi et al. / J. Nig. Soc. Phys. Sci. 4 (2022) 193–204 203 cells of humans under the influence of snake antivenom. In the future, we intend to use delay differential equations to model the impact of early and late treatment. Future studies may also be directed towards modeling the effects of other types of snake venom toxins, such as neurotoxins and anticoagulant toxins, on their target cells. Acknowledgments The authors are grateful to the anonymous reviewers and the handling editor for their constructive comments, which have undoubtedly improved the quality and clarity of the manuscript. We would also like to thank Dr. Abubakar S. Balla and Dr. Agom D. Ibrahim of Snakebite Treatment and Research Hospi- tal Kaltungo, Gombe state, for providing us with useful infor- mation on the management of snakebite envenoming patients. References [1] J. M. Gutierrez, R. D. G. Theakston & D. A. Warrell, “Confronting the neglected problem of snake bite envenoming: the need for a global part- nership”, PLoS Med 3 (2006) 6. [2] J. M. Gutierrez, D. Williams, H. W. Fan & D. A. Warrell, “Snakebite en- venoming from a global perspective: the need of an integrated approach”, Toxicon 56 (2010). [3] A. Kasturiratne, A. R. Wickremasinghe, N. de Silva, N. K. Gunawardena, A. Pathmeswaran, R. Premaratna, L. Savioli, D. G. Lallo & H. J. de Silva, “The global burden of snakebite: a literature analysis and modelling based on regional estimates of envenoming and deaths”, PLoS Med 5 (2008) 11. [4] J.P. Chippaux, “Snake-bites: appraisal of the global situation”, Bull World Health Organ 76 (1998) 5. [5] S. S. Williams, C. A. Wijesinghe, S. F. Jayamanne, N. A. Buckley, A. H. Dawson, D. G. Lalloo & H. J. de Silva, “Delayed psychological morbidity associated with snakebite envenoming”, PLoS Negl Trop Dis 5 (2011) 8, doi:10.1371/journal.pntd.0001255. [6] J. P. Chippaux,“Estimate of the burden of snakebites in sub-Saharan Africa: a meta-analytic approach”, Toxicon 57 (2011). [7] R. A. Harrison, A. Hargreaves, S. C. Wagstaff, B. Faragher & D. G.Lalloo, “Snake envenoming: a disease of poverty”, PLoS Negl. Trop. Dis. 3 (2009). [8] B. Mohapatra, D. A. Warrell, W. Suraweera, P. Bhatia, N. Dhingra, R. M. Jotkar, P. S. Rodriguez, K. Mishra, R. Whitaker, P. Jha & J. O. Gyapong, “Snakebite Mortality in India: A Nationally Representa- tive Mortality Survey”, PLoS Neglected Tropical Diseases 5 (2011) 4, doi:10.1371/journal.pntd.0001018. [9] Animal bites: Fact sheet N◦373, World Health Organization. https://www.who.int/mediacentre/factsheets/fs373/en/. Accessed: Febru- ary, 2021. [10] J. P. Chippaux, “Snakebite envenomation turns again into a neglected tropical disease”, J Venom Anim Toxins Incl Trop Dis. (2017). [11] J. M. Gutierrez, J. J. Calvete, A. G. Habib, R. A. Harrison, D. G. Williams & D. A. Warrell, “Snakebite envenoming”, Nat Rev Dis Primers 3 (2017) 17063. [12] World Health Organization, Prevalence of snakebite envenoming, https://www.who.int/snakebites/epidemiology/en/. Accessed: February, 2021. [13] World Health Organization, Animal bites, https:www.who.int/en/news- room/fact-sheets/detail/animal-bites. Accessed March, 2021. [14] World Health Organization, Regional Office for Africa. Guidelines for the prevention and clinical management of snakebite in Africa, (2010), http://apps.who.int/medicinedocs/ documents/s17810en/s17810en.pdf. Accessed: March, 2021. [15] B. S. Naik,“Dry bite in venomous snakes: A review” Toxicon, 133 (2017) 63, doi:10.1016/j.toxicon.2017.04.015. [16] B. A. Young & K. Zahn, “Venom flow in rattlesnakes: mechan- ics and metering”, Journal of Experimental Biology, 204 (2001) 24, doi:10.1242/jeb.204.24.4345. [17] D. A. Warrell, C. Bon & M. Goyffon, “Envenomings and Their Treat- ments. In: Clinical features of envenoming from snake bites”, Fondation Marcel Mérieux: Lyon, (1996). [18] L. A. Ribeiro, G. Puorto & M. T. Jorge, “Bites by the colubrid snake Philodryas olfersii. a clinical and epidemiological study of 43 cases”, Toxicon 37 (1999). [19] B. Alina, Live Science. Snake Facts and Types of Snakes. https://www.livescience.com/27845-snakes.html. Accessed March, 2021. [20] R. E. Hill & S. P. Mackessy, “Characterization of venom (Duvernoy’s secretion) from twelve species of colubrid snakes and partial sequence of four venom proteins”, Toxicon 38 (2000). [21] M. E.Peichoto, P. Teibler, R. Ruı́z, L. Leiva & O. Acosta, “Systemic pathological alterations caused by Philodryas patagoniensis colubrid snake venom in rats”. Toxicon 48 (2006). [22] S. Weinstein & D. Keyler, “Local envenoming by the Western hognose snake (Heterodon nasicus): A case report and review of medically signif- icant Heterodon bites”, Toxicon 54 (2009). [23] D. A. Warrell, N. M. D.Davidson & B. M. Greenwood, “Poisoning by bites of the saw-scaled or carpet viper (Echis carinatus) in Nigeria”, Q J Med. 46 (1977). [24] R. Bailey, How Does Snake Venom Work? https://www.thoughtco.com/how-snake-venom-works-4161270. Ac- cessed: February, 2021. [25] G. Leon, L. Sánchez, A. Hernández, M. Villalta, M. Herrera, A. Segura, R. Estrada & J. M. Gutiérrez, “Immune Response towards Snake Venoms. Inflammation and Allergy - Drug Targets”, Journal of Inflammation and Allergy Drug Target 10 (2011) 5. [26] B. S. Gold, C. D. Richard & A. B. Robert, “Bites of ven- omous snakes”, The New England Journal of Medicine, 5 (2002), doi:10.1056/NEJMra013477. [27] B. J. Daley & J. Torres, “Venomous snakebites”, A Journal of Emergency Medical Services 39 (2014) 6. [28] United State National Institute for Occupational Safety and Health, Venomous Snakes, Accessed: February, 2021, https://www.cdc.gov/niosh/topics/snakes/. [29] M. Peden & World Health Organization, “World report on child injury prevention”, World Health Organization, (2008), https://apps.who.int/iris/handle/10665/43851, Accessed May, 2021. [30] World Health Organization, Snakebite envenoming: a strategy for preven- tion and control, Geneva, (2019). [31] F. Chappuis, S. K. Sharma, N. Jha, L. Loutan & P. A. Bovier, “Protection against snake bites by sleeping under a bed net in southeastern Nepal”, Am. J. Trop. Med. Hyg. 77 (2007). [32] H. A. de Silva, A. Pathmeswaran, C. D. Ranasinha, S. Jayamanne, A. Hit- tharage, R. Kalupahana, G. A. Ratnatilaka, W. Uluwatthage, J. K. Aron- son, J. M. Armitage, D. G. Lalloo & H. J. de Silva, “Low-dose adrenaline, promethazine, and hydrocortisone in the prevention of acute adverse re- actions to antivenom following snakebite: a randomised, double-blind, placebo-controlled trial”, PLoS Med, 8 (2011) 5, doi: 10.1371/journal. [33] J. M. Gutierrez, L. Bruno, L. Guillermo, R. Alexandra, C. Fernando & A. Yamileth, “Trends in Snakebite Envenomation Therapy: Scientific, Tech- nological and Public Health Considerations. Current Pharmaceutical De- sign”, 13 (2007) 28, doi:10.2174/138161207782023784. [34] A. G. Habib, U. I. Gebi & G. C. Onyemelukwe, “Snake bite in Nigeria”, Afr J Med Med Sci. 30 (2001) 3. [35] A. G. Habib, “Public health aspects of snakebite care in West Africa: per- spectives from Nigeria”, J Venom Anim Toxins Incl Trop Dis, 19 (2013) 27. [36] A. G. Habib,“Effect of pre-medication on early adverse reactions follow- ing antivenom use in snakebite: a systematic review and meta-analysis”, Drug Saf. 34 (2011)10, doi: 10.2165/11592050-000000000-00000. [37] K. Maduwage, “Snakebite coagulopathy: controversies in understanding and management”, Sri Lanka Journal of Medicine 26 (2017) 2. [38] P. Malasit, D. A. Warrell, P. Chanthavanich, C. Viravan, J. Mongkolsapaya, B. Singhthong & C. Supich, “Prediction, preven- tion, and mechanism of early (anaphylactic) antivenom reactions in victims of snake bites”, British Medical Journal, 292 (1986), DOI: 10.1136/bmj.292.6512.17 [39] V. M. Morais & H. Massaldi,“Snake antivenoms: adverse reactions and production technology”,J Venom Anim Toxins incl Trop Dis. 15 (2009) 203 Abdullahi et al. / J. Nig. Soc. Phys. Sci. 4 (2022) 193–204 204 1. [40] Neglected tropical diseases: Snakebite. World Health Organization, https://www.who.int/health-topics/snakebite. Accessed: March, 2021. [41] United State National Institute for Occupational Safety and Health, “Ven- omous Snakes”, https://www.cdc.gov/niosh/index.htm, Accessed: April, 2021. [42] A. G. Habib, A. Kuznik, M. Hamza, M. I. Abdullahi, B. A. Chedi, J-P. Chippaux, & D. A. Warrell, “Snakebite is under appreciated:appraisal of burden from West Africa”, PLoS Negl.Trop. Dis. 9 (2015). [43] H. A. Reid & R. D. G. Theakston, “The Management of Snake bite”, Bulletin of the World Health Organization 61 (1983) 6. [44] M. A. Syed, A. A. Mohib, M. Jyotsna, C. Adarash & P. Jyotishka, “Emer- gency treatment of a snake bite: Pearls from literature”, J Emerg Trauma Shock 1 (2008) 2. [45] S. R. Vijeth, T. K. Dutta, J. Shahapurkar & A. Sahai, “Dose and frequency of anti-snake venom injection in treatment of Echis carinatus (saw-scaled viper) bite”, J Assoc Physicians India 48 (2000) 2. [46] D. A. Warrell, Clinical toxicology of snakebite in Africa and the Mid- dle East/Arabian peninsula, Asia In Meier J,White J eds. Handbook of Clinical Toxicology of Animal Venoms and Poisons. Florida, CRC Press, (1995). [47] D. Williams, J. M. Gutierrez, R. Harrison, D. A. Warrell, J. White, K. D. Winkel & P. Gopalakrishnakone, “The Global Snake Bite Initiative: an antidote for snake bite”, Lancet (London, England) 375 (2010) 9708, doi.org/10.1016/s0140-6736(09)61159-4 [48] D. A. Warrell, “Snakebite”, Erratum in: Lancet 375 (2010) 9708, doi: 10.1016/S0140-6736(09)61754-2 [49] Z. U. Abdullahi, S. S. Musa, D. He & U. M. Bello, “Antiprotozoal Effect of Snake Venoms and their Fractions: A Systematic Review”, Pathogens 10 (2021) 1632, doi.org/10.3390/pathogens10121632 [50] T. A. Reeks, B. G. Fry & P. F. Alewood, “Privileged frameworks from snake venom. Review”, Cell. Mol. Life Sci. (2015). DOI 10.1007/s00018- 015-1844-z [51] S. A. Zainal Abidin, Y. Q. Lee, Y. O. Lekhsan & N. Rakesh, “ Malaysian Cobra Venom: A Potential Source of Anti-Cancer Therapeutic Agents”, Toxins, 11 (2019) c2. [52] A. T. Mohamed, S. A. Garcia & J. D. Stockand, “Snake Venoms in Drug Discovery: Valuable Therapeutic Tools for Life Saving”, Toxins 11 (2019) 10, doi.org/10.3390/toxins11100564 [53] B. Lomonte, J. Fernández, L. Sanz, Y. Angulo, M. Sasa, J. M. Gutiérrez & J. J. Calvete, “Venomous snakes of Costa Rica: biological and med- ical implications of their venom proteomic profiles analyzed through the strategy of snake venomics”, J Proteomics 13 (2014) 105, doi: 10.1016/j.jprot.2014.02.020. [54] T. T. Yusuf, A. Abidemi, A. S. Afolabi & E. J. Dansu, “Optimal Control of the Coronavirus Pandemic with Impacts of Implemented Control Mea- sures”, Journal of the Nigerian Society of Physical Sciences 4 (2022) 1, doi: 10.46481/jnsps.2022.414. [55] A. S. Hassan & N. Hussaini, “Analysis of an HIV - HCV simultaneous infection model with time delay”, Journal of the Nigerian Society of Phys- ical Sciences 3 (2021) 1, doi.org/10.46481/jnsps.2021.109 [56] S. A. Ayuba, I. Akeyede, & A. Olagunju, “Stability and Sensitivity Analysis of Dengue-Malaria Co-Infection Model in Endemic Stage”, Journal of the Nigerian Society of Physical Sciences 3 (2021) 2, doi.org/10.46481/jnsps.2021.196 [57] I. Babaji Muhammad & S. Usaini, “Dynamics of Toxoplasmosis Disease in Cats population with vaccination”, Journal of the Nigerian Society of Physical Sciences 3 (2021) 1. [58] S. A. Abdullahi, A. G. Habib, & N. Hussaini, “Control of snakebite en- venoming: A mathematical modeling study”, PLOS Neglected Tropical Diseases 15 (2021) 8, doi.org/10.1371/journal.pntd.0009711(2021). [59] F. Y. Eguda, A. James & S. Babuba. “The Solution of a Mathematical Model for Dengue Fever Transmission Using Differential Transformation Method”, Journal of the Nigerian Society of Physical Sciences 1 (2019) 3, doi:10.46481/jnsps.2019.18. [60] L. Adamu & N. Hussaini, “An Epidemic Model of Zoonotic Visceral Leishmaniasis with Time Delay”, Journal of the Nigerian Society of Phys- ical Sciences 1 (2019) 1, doi: 10.46481/jnsps.2019.5. [61] A. Goyal, L. E. Laura & A. S. Perelson, “Within-host mathematical mod- els of hepatitis B virus infection: Past, present, and future. Current Opin- ion in Systems Biology”, (2019), doi.org/10.1016/j.coisb.2019.10.003. [62] H. Dahari, J. E. Layden-Almer, E. Kallwitz, R.M. Ribeiro, S. J. Cotler, T. J. Layden & A. S. Perelson, “ A mathematical model of hepatitis C virus dynamics in patients with high baseline viral loads or advanced liver disease”, Gastroenterology, 136 (2009) 4, doi.org/10.1053/j.gastro.2008.12.060 [63] O. Sharomi & A. B. Gumel, “Mathematical study of in-host dynamics of Chlamydia trachomatis”, IMA Journl of Applied Mathemtics 77 (2012) 2, doi:10.1093/imamat/hxq057 [64] L. M. de Freitas, T. U. Maioli, H. A. L. de Ribeiro, P. Tieri & F. Castiglione,“A mathematical model of Chagas disease in- fection predicts inhibition of the immune system”, [IEEE 2018 IEEE International Conference on Bioinformatics and Biomedicine (BIBM) - Madrid, Spain (2018.12.3-2018.12.6)] 2018 IEEE Interna- tional Conference on Bioinformatics and Biomedicine (BIBM), (2018), doi:10.1109/BIBM.2018.8621389. [65] A. Ishaiku & N. Hussaini, “Modelling Effects of HCV Plus-Strand RNA Influx into a Cell during HCV replication”, J.Indones. math. Soc. 26 (2020) 01. [66] P. Ngina, R. W. Mbogo & L. S. Luboobi, “The In-Vivo Dynamics of HIV Infection with the Influence of Cytotoxic T Lymphocyte Cells”, Interna- tional Scholarly Research Notices (2017), doi:10.1155/2017/2124789 [67] R. M. Kini & H. J. Evans, “A model to explain the pharmacological effects of snake venom phospholipases A2”, Toxicon (1989) 27. [68] S. Stagg, Pharmacokinetics of Snake Bites, http://www.umich.edu/ elements/fogler&gurmen/html/web−mod/cobra/avenom.htm. Accessed: January, 2019. [69] P. P. Tanos, G. K. Isbister, D. G. Lalloo, C. M. Kirkpatrick & S. B. Duffull, “A model for venom-induced consumptive coagulopathy in snake bite”, Toxicon 52 (2008). [70] The Lancet, “The Snakebite envenoming: a priority neglected tropical dis- ease”, Lancet, 1 (2017) 2, doi: 10.1016/S0140-6736(17)31751-8. PMID: 28677550. [71] J. Horky, J. Vacha & V. Znojil, “Comparison of life span of erythrocytes in some inbred strains of mouse using 14C-labelled glycine”, Physiol. Bohemoslovaca 27 (1977). [72] J. W. Goodman & L. H. Smith, “Erythrocyte life span in normal mice and in radiation bone marrow chimeras”. Am. J. Physiol. 200 (1961). [73] V. Lakshmikanthan, S. Leela & A. A. Martyniuk, A Stability Analysis of Nonlinear Systems, CRC Press, (1989). [74] K. Okuneye & A. B. Gumel, “Analysis of a temperature and rainfall- dependent model for malaria transmission dynamics”, Mathematical Bio- sciences 000 (2016). [75] Y. Dumont & F.Chiroleu, “Vector control for the Chikungunya disease”, Mathematical Biosciences and Engineering 7 (2010). [76] S. M. Sabiu & N. Hussaini, “Modeling the Dynamics of Dengue Virus in Human and Mosquito Populations”, Journal of Nigerian Association of Mathematical Physics 43 (2017). [77] J. P. LaSalle, The Stability of Dynamical Systems, Regional Conference Series in Applied Mathematics, SIAM Philadephia, (1976). [78] M. Lebois & E. C. Josefsson, “Regulation of Platelet lifes- pan by apoptosis. Platelets”, Journal of Platelets Sep (2016) 6, doi:10.3109/09537104.2016.1161739. [79] H. S. Bawaskar, “Snake venoms and antivenoms: Critical supply issues”, J Assoc Phys India 52 (2004). [80] A. G. Habib, M. Lamorde, M. M. Dalhat, Z. G. Habib & A. Kuznik, “Cost-effectiveness of Antivenoms for Snakebite Envenoming in Nigeria”, PLoS Negl Trop Dis 9 (2015) 1, doi.org/10.1371/journal.pntd.0003381 204