Fractional-Order Modelling and Optimal Control of Cholera Transmission

A Caputo-type fractional-order mathematical model for"metapopulation cholera transmission"was recently proposed in [Chaos Solitons Fractals 117 (2018), 37--49]. A sensitivity analysis of that model is done here to show the accuracy relevance of parameter estimation. Then, a fractional optimal control (FOC) problem is formulated and numerically solved. A cost-effectiveness analysis is performed to assess the relevance of studied control measures. Moreover, such analysis allows us to assess the cost and effectiveness of the control measures during intervention. We conclude that the FOC system is more effective only in part of the time interval. For this reason, we propose a system where the derivative order varies along the time interval, being fractional or classical when more advantageous. Such variable-order fractional model, that we call a 'FractInt' system, shows to be the most effective in the control of the disease.


Introduction
Fractional calculus is an old subject that raised as a consequence of a pertinent question that L'Hôpital asked Leibniz in a letter about the possible meaning of a derivative of order 1/2. Recently, many researchers have focused their attention in modelling real-world phenomena using fractional-order derivatives. The dynamics of those problems have been modelled and studied by using the concept of fractional-order derivatives. Such problems appear in, for example, biology, physics, ecology, engineering, and various other fields of applied sciences, see, e.g., [1][2][3].
Cholera is a gastroenteritis infection, contracted after consuming an infectious dose or inoculum size of the pathogenic vibrio cholerae [4]. The mode of transmission consists of two pathways: the primary route, where individuals consumes the pathogen from vibrio contaminated water and seafood; the secondary route being characterised by individuals consuming unhygienic or soiled food that is infested with pathogenic vibrios from an infected person. This secondary route of transmission is also commonly referred to as person-to-person contact [4]. Cholera infection has affected many parts of the world. However, its devastating force has been more pronounced in impoverished communities [5,6].
Cholera is one of the most studied infections in recent years. Mathematical models are used to study and understand the dynamics of infection, as well as offer suggestions toward its control, see, e.g., [4,[7][8][9][10][11][12]. Most models used in the study of Cholera have been based on systems of integer-order ordinary differential equations. However, those models do not fully account for memory, as well as non-local properties exhibited by the epidemic system.

arXiv:2112.03813v1 [math.OC] 7 Dec 2021
Non-local behaviour asserts that the subsequent state of the model depends on both the current and historical states. Fractional order differential systems have been proposed as more suitable to describe epidemic dynamics of diseases [13][14][15].
Recently, a fractional order differential system was used in the study of cholera transmission in adjacent communities [10]. Individuals in adjacent communities often frequent their home range, an aspect which constitutes memory that is present in such a model. Here, we start to do a sensitivity analysis to the fractional-order model of [10] in order to determine which model parameters are most influential on the disease dynamics. After that, fractional optimal control (FOC) is applied, as a generalisation of the classical optimal control system of [11]. In contrast with other works, for example [16], we could not find a derivative order for which the FOC system is more effective. However, we noticed that the FOC system is more effective in part of the time interval. Hence, we propose a system where the derivative order varies along the interval, being fractional or classical when more advantageous. Such system, a variable-order fractional one, is named here a FractInt system, shown to be useful in the control of the disease.
The paper is organised as follows. In Section 2, the fractional order model formulation is presented. Our main results are then given in Section 3: sensitivity analysis of the parameters of the model, taking into account the derivative order (Section 3.1); fractional optimal control of the model (Section 3.2); numerical simulations and cost-effectiveness of the fractional model (Section 3.3); and the new variable-order FractInt system (Section 3.3.2). We end with Section 4 of conclusions.

Fractional-Order Cholera Model
A metapopulation model for cholera transmission is considered, dividing the population into mutually exclusive distinct groups and using deterministic continuous transitions between those groups, also known as states. The model describes the dynamics of a population exposed to infection by the pathogen vibrio cholerae. The human population is divided into three compartments: susceptible individuals (S); infectious individuals (I); and recovered individuals (R). The Caputo fractional-order system of differential equations is as follows [10]: for the first sub-population, and 3 of 17 for the second sub-population, where C 0 D α t denotes the left Caputo fractional order derivative of order α [1], 0 < α 1, We note that the equations of model (1)-(2) do not have appropriate time dimensions. Indeed, on the left-hand side the dimension is (time) −α while on the right-hand side the dimension is (time) −1 (see, e.g., [17,18] for more details). Therefore, we claim that the accurate way of writing system (1) is for the first sub-population, and the accurate way of writing system (2), with Q 1 = µ α

Main Results
We begin by doing a sensitivity analysis to the parameters of the model, in order to identify those for which a small perturbation leads to relevant quantitative changes in the dynamics.

Sensitivity Analysis
Two distinct ways to compute the basic reproduction numbers, R 01 and R 02 , of the two sub-populations of the model, are available in [10,11]. To know which one is proper, we determine them by using the next-generation matrix method [19]. We obtain the community specific reproduction numbers as for the first community, and for the second community, where .
The values of used parameters are presented in Table 1, which were taken from [10], with exception of π 1 , π 2 and 2 . The first two parameters are equal and are defined as which is bigger than the values proposed in [10], in order to ensure an endemic scenario (R 01 > 1 and R 02 > 1). Indeed, it is the presence of an endemic situation that motivates us, in Section 3.2, to apply optimal control theory to tackle the cholera problem. The latter parameter, We start by considering that rates u, v, and m are all zero in the sensitivity analysis, unless specified otherwise. Definition 1 (See [27,28]). The normalised forward sensitivity index of R 0i , which is differentiable with respect to a given parameter p, is defined by Table 2 presents the values of the sensitivity index of the parameters of the model, obtained by the normalised sensitivity index (7), for the classical case (α = 1) of connected communities. These values have a meaning. For instance, Υ R 01 1 = +0.999 means that increasing (decreasing) 1 by a given percentage increases (decreases) always R 01 by nearly that same percentage. Sensitive parameters should be carefully evaluated, once a small perturbation in such parameter leads to significant quantitative changes. On the other hand, the estimation of a parameter with a small value for the sensitivity index does not require as much attention to evaluate, because a small perturbation in that parameter leads to small adjustments [29]. According with Table 2, we should pay special attention to the estimation of sensitive parameters 1 and 2 . In contrast, the estimation of K, µ p , β 1 , β 2 , σ 1 , σ 2 , g 1 , and g 2 do not require as much attention because of its low sensitivity. The missing parameters are those whose index value is zero. Table 2. Sensitivity of R 01 (top) and R 02 (bottom), evaluated for the parameter values given in Table 1 with For the fractional model, the sensitivity index depends on the derivative order α. We can see this in Figure 2, where: (a) the impact of variation of α in the sensitivity index of β 1 is displayed for the first community; (b) the impact of variation of α in the sensitivity index of 2 is exhibited for the second community. The graphics of the other parameters are not shown because they exhibited similar behaviours. The parameters whose index value for α = 1 is close to zero, as β 1 , do not vary much if we consider lower values of α, as we see in Figure 2a. On the other hand, parameters whose index value in Table 2 are not as close to zero as the previous one, as 2 , vary significantly if we consider lower values of α, as we see in Figure 2b, and their sensitivity decreases with the decrease in α.  Figure 3 presents the evolution of the sensitivity index for the specific reproduction numbers of the two communities, R 01 and R 02 , with the variation of the derivative order. We see that the evolution of the sensitivity index is analogous for the two reproduction numbers.

Fractional Optimal Control of the Model
Modelling dynamic control systems optimally is a very important issue in applied sciences and engineering [30]. In this section, our aim is to minimise the number of cholera infected individuals and, simultaneously, to reduce the associated cost. This is achieved trough: (i) the use of vaccination into communities, as an effective time-dependent measurecontrol v(t); (ii) the use of clean treated water, a preventive measure control u(t); (iii) and the implementation of proper hygiene, another preventive measure control m(t), in order to control person-to-person contact. Thus, we consider the following fractional optimal control problem: subject to system (3)-(4) with given initial conditions Note that we are using a quadratic cost functional on the controls, as an approximation to the real non-linear functional. Indeed, as in [31], our optimal control problem depends on the assumption that the cost takes a non-linear form. The parameters 0 < k 1 , k 2 , k 3 , k 4 , k 5 < +∞ are positive weights and t f is the duration of the control program. In addition, k 3 , k 4 , and k 5 represent the costs of applying controls efforts u, v, and m, respectively. The set of admissible control functions is The Pontryagin maximum principle (PMP) for fractional optimal control is used to solve the problem [32,33]. The Hamiltonian of the resulting optimal control problem is defined as and the adjoint system asserts that the co-state variables ξ i (t), i = 1, . . . , 8, verify with t = t f − t. In turn, the optimality conditions of PMP establish that the optimal controls u * , v * , and m * are defined by u * = min max 0, In addition, the following transversality conditions hold:

Numerical Results and Cost-Effectiveness Analysis
We start by calculating the relevance of the three measures used in the control of the disease. We do it by using the sensitivity index, presented in Definition 1, as proposed in [34]. In this case, the sensitivity indices are presented as functions of the control parameters in Figure 4, using the parametric values from Table 1 and considering the classical model (that is, α = 1). The resulting graphics show that: (a) the curve of vaccination is the one that most rapidly moves away from zero, meaning that the vaccination program has a big impact for small rates of application; (b) proper hygiene measures are also important in the control of cholera, having a bigger impact for greater rates of application of it, being the only control that has precisely the same impact in both communities; (c) the domestic water treatment is useless in the control of cholera transmission, when used simultaneously with the two previous controls. So, in what follows, the third control, variable u, is ignored.  One can deal with many different problems arising in different fields of sciences and engineering by applying some appropriate discretisation [35]. Here, the Pontryagin Maximum Principle is used to numerically solve the optimal control problem, as discussed in Section 3.2, both in classical and fractional cases, using the predict-evaluate-correct-evaluate (PECE) method of Adams-Basforth-Moulton [36], coded by us in MATLAB. Firstly, we solve system (3)-(4) by the PECE procedure, with the initial values for the state variables as given in Table 3 and a guess for the controls over the time interval [0, t f ], that way obtaining the values of the state variables. Similarly to [16], a change of variable is applied to the adjoint system and to the transversality conditions, obtaining the fractional initial value problem (12)- (14). Such IVP is also solved with the PECE algorithm, and the values of the co-state variables ξ i , i = 1, . . . , 8, are determined. The controls are then updated by a convex combination of the controls of previous iteration and the current values, computed according to (13). This procedure is repeated iteratively until the values of all the variables and the values of the controls are almost identical to the ones of the previous iteration. The solutions of the classical model were successfully confirmed by a classical forward-backward scheme, also implemented by us in MATLAB.
In the numerical experiments, we consider weights k 1 = 4, k 2 = 2.4, k 3 = 1.6, and k 4 = k 5 = 1. These values have the same relations between them than the homonyms weights in [11], with the numerical advantage of being closer to the value one. We also use v max = m max = 1, while the other parameters are fixed according to Table 1.
According with Figure 1, R 01 > 1 and R 02 > 1 for α 0.68 (endemic scenario). In order to be able to compare the FOCP results for several derivative orders, we consider the initial conditions, given by Table 3, which correspond to the non-trivial endemic equilibrium for system (3)-(4) for the classical cholera model (α = 1). Table 3. Initial conditions for the fractional optimal control problem of Section 3.2 with parameters given by Table 1, corresponding to the endemic equilibrium of cholera model (3)-(4) with classical derivative order. Without loss of generality, we consider the fractional order derivatives α = 1.0, 0.9 and 0.8. In Figures 5-7, we find the solutions of the fractional optimal control problem for those values of α. We can see that a change in the value of α corresponds to significant variations of the state and control variables. Beyond those values of α, others values were also tested, but the results do not changed qualitatively.
The efficacy function [37], exhibited in Figure 8, is defined as where i * (t) = I * 1 (t) + I * 2 (t) is the optimal solution associated with the fractional optimal control problem and i(0) = I 1 (0) + I 2 (0) is the correspondent initial condition. This function measures the proportional variation in the number of infected individuals, of both communities, after the application of the control measures, {v * , m * }, by comparing the number of infectious individuals at time t with its initial value i(0). We observe that the graphic of F(t) exhibits the inverse tendency of infected individuals curves, growing and reaching the maximum at the end of the time interval.
To assess the cost and the effectiveness of the proposed fractional control measure during the intervention period, some summary measures are presented. The total cases averted by the intervention during the time period t f is defined in [37] by where i * (t) is the optimal solution associated with the fractional optimal controls and i(0) is the correspondent initial condition. Note that this initial condition is obtained as the equilibrium proportion i of systems (3)-(4), which is independent of time, so that t f i(0) = t f 0 i dt represents the total infectious cases over a given period of t f days. Effectiveness is defined as the proportion of cases averted on the total possible cases under no intervention [37]: The total cost associated with the intervention is defined in [37] by where s * (t) = S * 1 (t) + S * 2 (t) and C i corresponds to the per person unit cost of the two possible interventions: (i) vaccination at any time t of susceptible individuals (C 1 ); and (ii) infected individuals practising proper hygiene (C 2 ). Following [37,38], the average cost-effectiveness ratio is given by The cost-effectiveness measures are summarised in Table 4. The results show the effectiveness of controls to reduce cholera infectious individuals and the leadership in doing so by the classical model (α = 1). Table 4. Summary of cost-effectiveness measures for classical and fractional (0 < α < 1) cholera disease optimal control problems. Parameters according to Tables 1 and 3 Table 1 and fractional order derivatives α = 1.0, 0.9 and 0.8.

Optimal Control Strategies and Cost-Effectiveness Analysis
Now, we analyse the cost-effectiveness of alternative combinations of two possible control measures: • strategy A-implementing the vaccination control, v; • strategy B-implementing proper hygiene control, m; • strategy C-implementing both controls, v and m.
To analyse the cost-effectiveness of the three alternative strategies A, B, and C, we use the Incremental Cost-Effectiveness Ratio (ICER) [38]. This ratio is used to compare the differences between costs and health outcomes of two alternative intervention strategies that compete for the same resources, being often described as the additional cost per additional health outcome. We start by ranking the strategies in order of increasing effectiveness, assessed by the total averted cases AV, defined in (16).
The ICER for the classical model (α = 1), is calculated as follows: Results are shown in Table 5. Strategy A is the most costly one, so we exclude this strategy from the set of alternatives. We align the remaining strategies by increasing effectiveness (AV) and recalculate the ICER: ICER(B) = ACER(B) = 0.216276 and ICER(C) = 3.210922. Hence, we conclude that strategy B (implementing only the control of hygiene measure, m) is the most cost-effective strategy. Table 5. Incremental cost-effectiveness ratio for strategies A, B, and C. Parameters according to Tables 1  and 3 with C 1 = C 2 = 1 and α = 1. The ICER was also computed for the fractional model, considering the same strategies, with the three derivative orders used previously: α = 0.9, α = 0.8 and α = 0.68. Even in those cases, we obtain the same conclusion: strategy B is the most cost-effective.

Strategies
Comparing the cost-effectiveness of strategy B for above derivative orders, we start by getting the values of ICER presented in Table 6. Once the values of Total Cost (TC) diminish with the decrease in the derivative order, proceeding as above, when the values of ICER were computed for the classical α = 1 model, we eliminate successively the scheme with the highest TC. Consequently, the "strategies" α = 1.0, α = 0.9, and α = 0.8 are excluded, by this order.
Our conclusion is: the most cost-effective scheme is strategy B with α = 0.68. The solution of the Fractional Optimal Control Problem, exhibited in Section 3.3 for some derivative order values, evidences that the fractional-order model can be more effective in part of the time interval-as shown in Figure 8-but classical model is more effective if we consider all the interval, as shown in Table 4. In this section, we consider that the derivative order varies along the interval, being fractional or classical when more advantageous. According with Figure 8, α should be fractional at beginning and become classical/one after a certain time. Such class of Variable-Order Fractional (VOF) systems [39] are baptised here as FractInt. The derivative order of the FractInt system varies according with where 0 < α 0 < 1. In practice, we noticed that the resulting system is as more effective as smaller the value of α 0 . In view of an endemic scenario, we consider in our simulations α 0 = 0.68, which the lowest value that can guarantee it. With respect to the switching time, the value considered is t = 7, which is the one to which it corresponds the maximum value of efficiency.
The FractInt system is numerically solved with the procedure described above, at the beginning of Section 3.3. In each iteration of the procedure, the predict-evaluate-correctevaluate method is applied successively to each one of the two initial value problemsassociated to the state system and to the adjoint system of the FOCP (8)-(13)-to perform the integration of them in two steps. These steps correspond to the two branches of the derivative order function, defined by (20). In such process, the final solution of first step (first branch) is the initial solution of the second step.
Solutions of the FractInt system are reproduced in Figures 9 and 10 along with solutions of the classical and fractional models (α = 1 and α = 0.68, respectively). We can see that the solutions of the FracInt system start to follow the solutions of the fractional model and end by following the solutions of the classical model. This behaviour can be observed in Figure 11, where the efficacy function for these three models is displayed.
The cost-effectiveness measures for the FractInt system are summarised in Table 7. Table 7. Cost-effectiveness measures for the FractInt system. Parameters according to Tables 1 and 3 with C 1 = C 2 = 1. Our results show that it is effective to control cholera infection through optimal control and that the FracInt model is more efficient than the classical one (cf. Table 4 and Figure 11), being the most effective model.

Conclusions
Nowadays, Cholera infection is still an healthcare problem in many parts of the world, namely in impoverished regions where it has devastating effects. Its control is the goal of many studies in last years. In this work, a fractional-order mathematical model for Cholera with two connected communities, proposed in [10], is studied and generalised. The two communities in appreciation are distinct and therefore they could not be molten in one unique community.
A sensitivity analysis of the model is done to show the importance of estimation of parameters. A fractional optimal control (FOC) problem is then formulated and numerically solved.
The relevance of studied controls is assessed using a cost-effectiveness analysis. Such analysis allow us to neglect the control u (domestic water treatment) since it proved to be useless when used in combination with remaining controls.
The numerical results show that the FOC system is more effective only in part of the time interval. Therefore, we propose a system where the derivative order varies along the time interval, being fractional or integer when more advantageous in the control of infection. Such variable-order fractional model, baptised here as FractInt, shows to be the most effective in the control of the disease.