Modelling the Role of Human Behaviour in Ebola Virus Disease (EVD) Transmission Dynamics

The role of human behaviour in the dynamics of infectious diseases cannot be underestimated. A clear understanding of how human behaviour influences the spread of infectious diseases is critical in establishing and designing control measures. To study the role that human behaviour plays in Ebola disease dynamics, in this paper, we design an Ebola virus disease model with disease transmission dynamics based on a new exponential nonlinear incidence function. This new incidence function that captures the reduction in disease transmission due to human behaviour innovatively considers the efficacy and the speed of behaviour change. The model's steady states are determined and suitable Lyapunov functions are built. The proofs of the global stability of equilibrium points are presented. To demonstrate the utility of the model, we fit the model to Ebola virus disease data from Liberia and Sierra Leone. The results which are comparable to existing findings from the outbreak of 2014 − 2016 show a better fit when the efficacy and the speed of behaviour change are higher. A rapid and efficacious behaviour change as a control measure to rapidly control an Ebola virus disease epidemic is advocated. Consequently, this model has implications for the management and control of future Ebola virus disease outbreaks.


Introduction
Population growth and density in health risk areas of our societies do not facilitate the management of disasters whose number keeps growing as time evolves [1]. When a catastrophe arises, people can exhibit controlled behaviour which can take the form of intelligent and reasoned reactions [1]. However, panic or lack of knowledge can lead to a less controlled behaviour demonstrated by consideration or automated behaviour for example [1]. The efficacy of behaviour change can then be evaluated by the level of control in the behaviour of the affected population. During EVD outbreak of 2014 − 2016, the poor economic situation of the affected countries and some cultural beliefs impacted the level of control in the populations' behaviour, especially at the early stage of the epidemic.
High population mobility across porous borders, cultural beliefs, community resistance, strikes by health per-sonnel, and messages fuelling despair are all behaviourrelated factors listed by the World Health Organization (WHO) which contributed to the rapid and invisible spread of EVD [2]. Culturally meaningful behaviours such as washing and touching the deceased and bush-meat consumption practised in West Africa, unfortunately, helped to fuel EVD spread [3]. The disease was new to the unprepared health system of the affected countries which detected it three months after its onset in Guinea. So, a faster and more efficient behaviour change of the health authorities and the population would have greatly limited EVD spread and burden in West Africa [4,5]. As a consequence, any EVD modelling framework that aims at controlling the disease must incorporate the efficacy and the speed of behaviour change and this is done in this manuscript. Mathematical models of Ebola virus disease (EVD) can be divided into two groups: models that do not account for infection due to EVD deceased and models that do. Models ignoring the deceased are classified as SEIR or SIR models; see for instance [6][7][8][9][10][11][12]. The incidence rate of these models is generally bilinear, i.e., of the form βSI. The formulation of such an incidence function increases the ease of analysing the models, including the global stability analysis of the fixed points using LaSalle's invariance principle [6,10,12]. However, neglecting the contribution of the deceased in the force of infection of EVD underestimates the number of Ebola cases [13].
Several models of EVD including a compartment for the deceased have been formulated and analysed; see for instance [14][15][16][17][18][19]. The mathematical analysis of the equilibria of such models mostly focuses on the determination of the reproduction number and local stability analysis of equilibria [17,18]. The proof of the global stability of the equilibria is often neglected because of the nonlinearity of the incidence rate. Rivers et al. [15] modelled the impact of interventions on an epidemic of Ebola using a model that accounts for infection during funerals (F). The total population N, in this case, was made up of susceptibles ðSÞ, exposed ðEÞ, infected ðIÞ, hospitalized ðHÞ, and recovered ðRÞ individuals. The incidence rate of the model which was given by ðβSI + β H SH + β F SFÞ/N, accounted for infection during hospitalisation and burials. A similar structure of the incidence rate was used by Djiomba and Nyabadza [17] to investigate optimal control of EVD with educational campaigns, active case finding, and pharmaceutical interventions and Xia et al. [20] who modelled the transmission dynamics of Ebola in Liberia.
More complex incidence rates or forces of infection are found in disease models, especially those that incorporate behaviour change. Stopping a disease outbreak, in general, relies on the technical implementation of control strategies like vaccination, treatment, or educational campaigns. In the absence of treatment, a change in human behaviour will be the only hope to stop the disease spread. So, human behaviour remains the mainstay of any disease eradication, including EVD [21].
Psychology and anthropology are not the only domains interested in human behaviour. Mathematical modelling has also focused on the effects of human behaviour on various diseases' evolution; see for instance [21][22][23][24][25]. However, the quantification of human behaviour remains a challenge, and its effects on disease evolution are often represented by probabilities or by a modified incidence rate of the disease; see for instance [22]. Del Valle et al. [26] formulated a model to investigate the effects of behavioural changes in a smallpox attack. They assumed that individuals within a community change their behaviour and move from a normal active group to a less active group, reducing their average number of contacts, in response to a high prevalence of smallpox in the community. The per-capita transfer rate between compartments was modelled by the function ϕ i = a i ðI n + I l Þ/ð1 + b i ðI n + I l ÞÞ per day with i = S, E, I and where n = normally active, l = lessa ctive, and a i and b i being positive constants that modulate the rate of change. Behavioural change was thus taken to vary with the ϕ i values and the transmission rate of smallpox. The force of infection for the normally active and less active was given by where A c and W were, respectively, the confined and quarantined individuals, j = ðn, lÞ. They concluded that, besides the standard intervention procedure that would affect the extent and duration of a smallpox epidemic, the reduction in contacts of people in response to information about the epidemic was another contribution to smallpox control [26]. Wang X et al. [27] modelled the influence of human behaviour on cholera dynamics with a disease prevalence contact rate β i = a i − b i m i ðIÞ, where a i is the contact rate without considering the influence of human behaviour, b i is the maximum reduced contact rate due to human behaviour, and m i ðIÞ = I/I + K i , i = 1, 2, 3 where K i is the halfsaturation constant function. The contact rate was driven by a saturation function of the Michaelis Menten type. They concluded that human behaviour can reduce the endemic and epidemic levels, cholera spreading speeds, and the risk of infection [27]. Fractional derivatives are also recently used in disease modelling [28,29]. A fractional model with behaviour is presented in [30]. This model uses fractional operators to investigate the asymptomatic behaviour of immunogenic tumour dynamics. The predictor-corrector method is adapted to the tumour-cell population dynamics and used to develop a tracking control method that helps to limit the growth of the tumour-cell population [30].
The impact of behaviour change on the spread of EVD with a linear force of infection has been done in [31]. Most of the EVD models that have been formulated to date are of the types SEIR, SIR, SIRF, SEIRF, and SEIHRF: In this paper, we consider a simple SIRD model, where D is the class of the deceased with a new nonlinear incidence function. In general, people abandon good behaviour when the perceived risk of a disease is reduced or because of limited funds to continue [32]. In the case of EVD, we assume that when the number of infected cases is well pronounced, the fear of contracting the disease obliges people to change their behaviour and a decline in the number of cases follows. Disease incidence or prevalence in a population has already been mentioned as a motivator of behaviour change during an epidemic for the case of influenza [26,32]. Exponential functions have been suggested to represent the probability of disease transmission in the case of alertness to disease by Lio et al. [33]. Fast et al. [34] used a decreasing exposition probability, due to behaviour change between susceptible and infected individuals, to assess the role of social mobilisation in EVD control in Lofa.
The absence of treatment for EVD is a source of fear within the affected population and motivates a change in their behaviour leading to reduced disease transmission, especially when the number of infected cases reaches unexpected high values. In this model, we innovatively introduce a nonlinear force of infection that considers the efficacy and human behaviour change which were not considered in previous models [4,5,34]. EVD started in a rural area where it had been unnoticed for a few months, causing the force of infection to initially grow exponentially and later declined as the epidemic progressed [35]. Concave exponential 2 Computational and Mathematical Methods in Medicine functions are thus suitable to represent the force of infection. Among all the existing exponential functions, we propose a new Ricker-type function to represent the force of infection of our model. The function innovatively takes into consideration the efficacy and speed of behaviour change. It models the infection growth in the presence of behaviour change and the model fits better to EVD data from the two countries considered; see also [36]. The effective transmission rate of EVD is represented by the parameter β and the efficacy of behaviour change is represented by a parameter p ðp ∈ ½0, 1Þ such that p = 0 corresponds to an uncontrolled behaviour that fuels EVD spread and p = 1 corresponds to a perfectly controlled behaviour that stops EVD spread. We also introduce a parameter K that represents the speed at which individuals change their behaviour. The parameter K helps to track how fast individuals change their behaviour in response to EVD. K = 0 corresponds to no behaviour change and greater values of K indicate a positive change in the affected population's behaviour. Total eradication of EVD corresponds to a zero force of infection for very high values of KðK ⟶ ∞Þ. We propose a force of infection given by where ε indicates the infectivity of the deceased when compared to that of the infected individuals, with ε > 1 and lim K⟶∞ λ = 0. This paper aims to study the dynamics of an EVD model with a nonlinear force of infection innovatively described by a Ricker function. We built the model with the assumption that behaviour change is motivated by disease incidence. We present the stability analysis of the model's steady states, and numerical simulations are used to confirm the modelling assumptions and to quantify the parameters modelling the efficacy and speed of behaviour change.
The paper is arranged as follows: in Section 2, we formulate the model. Section 3 is reserved for the model analysis, and in Section 4, we focus on the steady states and their global properties. Section 5 is reserved for the numerical simulations and the last section contains some concluding remarks.

Model Formulation and Equations
A deterministic model, consisting of individuals with different EVD statuses, is formulated to represent the population dynamics when there is a change of behaviour due to the high prevalence of the disease. Susceptible individuals (in compartment S) represent individuals that are at risk of contracting the infection through contact with the infected (in compartment I) and the deceased (in compartment D). Susceptible individuals are recruited into a heterogeneous population at a constant rate Λ. Infected individuals may recover at the rate α or may succumb to the Ebola virus and die at a rate φ. Individuals in each compartment are assumed to die naturally at a rate μ, and dead bodies of EVD-infected individuals deceased are disposed at a rate ρ. Dead bodies of infected individuals are assumed to be infectious. We thus assume here that infected individuals who do not die from EVD contribute to EVD spread and we set δ = φ + μ [37]. The total population size is given by Figure 1 represents the movements between different classes, as their infection status with respect to the disease changes.

Model Properties
3.1. Existence and Uniqueness of Solutions. The right-hand side of the system of differential equations (4)- (7) is made of Lipschitz continuous functions. According to Picard's existence theorem, with given initial conditions, the solutions of our system exist and they are unique [38].

Theorem 1. The system makes biological sense in the region
which is attracting and positively invariant for the flow of system (4)- (7).
Integrating (7) gives DðtÞ decreases to DðtÞ ≥ Λ δ/ρ μ when t ⟶ ∞ and enters Ω. Besides, any sum or difference of variables in Ω with positive initial values will remain in Ω or in a neighbourhood of Ω. Thus, Ω is bounded, positively invariant, and attracting for the flow of system (4)-(7).

Positivity of Solutions
Theorem 2. The existing solutions of our system (4)-(7) are all nonnegative.
Equations (11)-(13) admit a disease free equilibrium (DFE) denoted by E 0 = ðΛ/μ, 0, 0Þ and the next generation matrix method (see [40]) yields a reproduction number R 0 given by R 0 = β ð1 − pÞΛ/μð1/Q 1 + ε ðδ/ρ Q 1 ÞÞ. From the expression of the reproduction number, we can state that a more effective behaviour change reduces the number of secondary EVD infections.
The endemic equilibrium E * = ðS * , I * , D * Þ is given by I * is solution of the equation obtained by substituting equations (14) and (15) into equation (12). We set MðI * Þ = β ð1 − pÞ S * I * ð1 + δ ε/ρÞ exp ½−I * ð1 + δ ε/ρÞK − Q 1 I * and find the roots in Θ of the function M. However, the exponential function does not facilitate the algebraic resolution of equation (16). We first rewrite the function M as MðI * Þ = ðΓðI * Þ − Q 1 Þ I * where ΓðI * Þ = β ð1 − pÞ S * ð1 + ε ðδ/ρÞÞ exp ½−I * ð1 + δ ε/ρÞK and I * = 0 is a root of M which corresponds to the DFE. The endemic equilibrium point (I * ) is the solution of ΓðI * Þ = Q 1 . We have and Γð0Þ = R 0 Q 1 and lim I * ⟶∞ ΓðI * Þ = 0. So ΓðI * Þ is a monotone decreasing and positive function which intercepts the line y = Q 1 only once when R 0 > 1. In fact, the graph of Γ can intercept the line y = Q 1 if Γð0Þ ≥ y ð0Þ because Γ is a decreasing function. Since Γð0Þ = R 0 Q 1 , Γð0Þ is above Q 1 when R 0 > 1. If R 0 < 1, Γð0Þ is below Q 1 , implying that we have no endemic equilibrium. We can then conclude that our model has a unique endemic equilibrium point when R 0 > 1.
We use the parameter values obtained by the fitting process in Figure 2(a) to illustrate the existence of a unique endemic equilibrium point in Figure 3.
where A and G are positive constants to be calculated with FðE * Þ = 0. The right-hand side of system (11)-(13) at equilibrium yields where g = δ/ρ. The derivatives of F with respect to each state variable are

Computational and Mathematical Methods in Medicine
The derivative of F with respect to time t, denoted by _ F, is Considering the system of equations (11)-(13) together with inequality (??) and the fact that E * is the global minimum of F, we obtain We set and rewrite (25) as where The expression of L is obtained by replacing Λ, Q 1 , and g in equation (25) by their expressions from equations (21).
Expanding the expression of Lðx, y, uÞ from system (28) and grouping the coefficients with the same variable give Since we already have −μ ððS − S * Þ 2 /SÞ ≤ 0 from the expression of _ F, it remains to prove that L ≤ 0 in order to get _ F ≤ 0. From the expression of L in (29), we will first set the terms containing variables and with nonnegative coefficients to zero in order to get rid of the positive and nonconstant part of L. The coefficients of y, u, x y, and u x are thus set to zero and solved for A and G. We obtain So L is negative and will be equal to zero if x = y = u = 1 which is in the set LaSalle's extension implies that each solution which intersects ℝ 3 + limits to the endemic equilibrium and E * is globally asymptotically stable on Θ (see [41]).
Global stability of the EE indicates that EVD persists as long as the value of the reproduction number is greater than one. Programs aiming at eradicating EVD should limit the disease transmission in such a way that the value of the reproduction number remains less than one.

Model Validation.
We fit the number of new EVD cases IðtÞ from our model to data from Liberia and Sierra Leone [42]. Table 1 gives the number of new Ebola cases recorded in two countries. The values of the parameters estimated 6 Computational and Mathematical Methods in Medicine through the least-squares fitting process are given in Table 2.
A comprehensive description of the least-squares fitting method is given in [19]. Figures 2(a) and 4(a) show that our model closely fits to the data and indicate that behaviour change during the EVD outbreak of 2014 − 2016 was motivated by an increased number of EVD cases. They also show that an efficacious behaviour change was necessary to stop EVD spread. Figures 2(b) and 4(b) show a poor fit of the model to data from Liberia and Sierra Leone when there is a slow and inef-ficacious behaviour change, resulting in an uncontrolled behaviour during EVD epidemic (K = 0:005, p = 0). The controlled or efficacious behaviour adopted by populations in order to limit EVD spread consisted of practising careful hygiene, avoiding contact with infected animals and humans, collaborating with case tracking, not fleeing from isolation areas, safely burying those who died from EVD, etc. [17,43]. Figures 2(b) and 4(b) predict an increased number of infected individuals (represented by the continuous line) when compared to Figures 2(a) and 4(a), thus    overestimating the number of EVD cases. This increase is due to a slower and nonefficacious behaviour change that caused the EVD epidemic to last longer and to claim more human lives. Behaviour change should be advocated at the beginning of an epidemic, to limit the death toll of the disease and its spread.
The values of the estimated parameters in Table 2 yield a reproduction number equal to 1:71 in Liberia and 2 in Sierra Leone when there is a rapid and controlled change in human behaviours. These values are comparable to those obtained by several authors. Rivers et al. [15] found that the value of the reproduction number for the two countries on an average is 2:2. Althaus [44] found R 0 = 1:59 for Liberia and R 0 = 2:53 for Sierra Leone in 2014. Xia et al. [20] found R 0 = 2:02 for Liberia when investigating the different transmission routes of EVD. Chowell and Nishiura [11] estimated the value of the reproduction number between 1 and 2 in Liberia and Sierra Leone, from March to August 2014. The difference between the estimated and cited values of the reproduction number might come from the fact that they are calculated during different periods of EVD epidemics and the models built by the authors have different structures and integrate different control measures in some instances.
Behaviour change influenced several factors such as the EVD transmission rate, death rate, recovery rate, and even the burial rate of the deceased. Avoiding contact between susceptible and infected individuals, for example, contributed to a reduction in the transmission rate of EVD. This explains why the value of β estimated in Table 2 when controlled behaviour is considered is less than the one from the fitting with slow behaviour change and uncontrolled behaviour. Some researchers concluded that the transmission rate of EVD is higher without behaviour change [15,20]. The recovery rate and death rate of EVD are influenced by the collaboration between populations and authorities through active case finding or hospitalisation for example. Recover-ing from EVD without treatment was very rare and can explain the low values of α estimated in Table 2. The case fatality rate of the last EVD epidemic was greater than 50% [15,20,37], indicating a high risk of death for EVD patients. This supports the estimated values of δ in Table 2 for both countries. In [45], the values of δ are greater without a controlled behaviour, highlighting the importance of adopting disease reducing behaviours like visiting the appropriate health centres to receive pharmaceutical interventions that help to increase the chances of recovering from EVD. Dead bodies of EVD deceased are twice more infectious than the living infected individuals according to the estimation of ε. Some authors found that this infectivity could be even four times higher for the deceased [20]. Quick and safe burials of EVD deceased have helped to limit EVD spread, and the high values of ρ from Table 1 indicate that a high burial rate helped to limit the disease spread. The estimated values of p indicate that a highly efficacious behaviour change, adapted to health authorities' instructions, for example, was necessary to halt EVD spread. The estimated values of K indicate how fast a behaviour change occurred during the last EVD epidemic in Liberia and Sierra Leone. However, these values are low and this implies that behaviour change was a slow process during the outbreak. This could explain the long duration of the outbreak and its extent in terms of the number of deaths. Cultural beliefs and poor economic conditions build mistrust between rural populations and authorities in West Africa [3]. This situation led to a slow application of instructions from authorities during the outbreak. So, building confidence and improving people's living conditions is essential in motivating people to quickly react in case of an epidemic of EVD.

Sensitivity
Analysis. The reproduction number R 0 is made of different parameters and their correlations to R 0 differ as well. Sensitivity analysis helps to assess the variability  [46]. We use Latin hypercube sampling (LHS) and partial correlation coefficient (PRCC) to explore the parameter space of the model with 1000 simulations [47]. We choose a uniform distribution to implement LHS sampling scheme. The PRCC value for a specific parameter is a Pearson correlation coefficient for the residuals from two regression models [46]. The most important parameters have a PRCC value greater than 0:5 or less than −0:5. Figure 5 represents the sensitivity analysis of R 0 . We observe in Figure 5 that the parameters with the most important positive correlations to R 0 are δ, Λ, and β and those with the negative correlations are ρ, μ, and p. The results show that when the transmission rate and death rate of EVD increase, the number of new infections increases since the deceased and infected individuals are both infectious. When the recruitment rate is increased, more individuals are exposed to the Ebola virus and this increases their chances of contracting the virus and contributes to the disease's spread. Natural death reduces the size of the population and thus limits the number of people who might be infected by EVD. Disposing of the corpses of EVD deceased rapidly and safely stops contamination during funerals and reduces the potential number of EVD cases. The negative correlation between R 0 and p is certainly because adopting a controlled behaviour that limits or stops EVD spread reduces the reproduction number and is represented by an increase of the value of p.

Conclusion
The global stability analysis of the steady states of a model of EVD with a nonlinear incidence function is presented in this paper. The global stability of equilibria is presented through suitably chosen Lyapunov functions. The nonlinear incidence function is chosen to represent the influence of human behaviour on EVD evolution. Data from Liberia and Sierra Leone are used in the fitting process whose results show that the model closely describes the evolution of EVD with human behaviour. Parameters obtained in the fitting process support the assumptions made in the development of the model that behaviour change is motivated by disease prevalence. As soon as an epidemic starts, individuals affected should trust the authorities and do as instructed to avoid deaths and disasters due to a slow reaction.
The work presented in this paper is not without shortcomings. The modelling of human behaviour through any mathematical function is a challenge as human behaviour is complex and unpredictable. A stochastic version of the incidence function may be worth looking at as an alternative albeit the challenges associated with its formulation and the mathematical analysis. The use of individual-based models or network models could potentially solve the inadequacies of dealing with deterministic models. We still however argue that this model presents some interesting mathematical results on the global stability of an EVD model with human behaviour and a nonlinear incidence rate.
The theory of fractional calculus is often used for modelling purposes in general [48][49][50]. In disease modelling, in particular, fractional derivatives which are used in modelling the dynamics of an infectious disease like COVID-19 [51] could also be applied to Ebola virus disease models in future works.

Disclosure
This paper is published in a dissertation [52].

Conflicts of Interest
The authors declare that they have no conflicts of interest.