COST-EFFECTIVENESS AND BACKWARD BIFURCATION ANALYSIS ON COVID-19 TRANSMISSION MODEL CONSIDERING DIRECT AND INDIRECT TRANSMISSION

Since it was first discovered in Wuhan, China, COVID-19 has continued to spread throughout the world. Since then, many research works have been conducted to understand the spread of COVID-19. In this article, we propose an epidemiological model to understand the spread of COVID-19, considering the saturated treatment rate, direct/indirect transmission, and optimal control problem to find the best strategy for the COVID-19 eradication program. The model constructed is based on a nonlinear system of ordinary differential equations. Analytical results regarding the basic reproduction number and all equilibrium points are obtained analytically. Our model shows a possibility of the existence of the COVID-19 endemic state such that even the basic reproduction number is less than unity. We also found that indirect transmission contributes to the increases in the basic reproduction number and also the occurrence of the multiple endemic states. An optimal control approach was applied to determine the best strategy for the COVID-19 eradication program. Three control parameters were considered in the model: medical mask, disinfectant, and medical treatment. A Pontryagin’s Maximum Principle was used to derive the optimal control characterization of the related model and was solved numerically using the forward-backward iterative method.


INTRODUCTION
At the end of 2019, the world was shocked by the discovery of a new variant of the coronavirus known as the COVID-19 [1]. The strongest suspicion is that the virus was first spread through animals and then transmitted to humans. Since it was first discovered in the animal market in Wuhan, China, the disease has continued to spread throughout the world. According to data from WHO until 15 July 2020, more than 13 million cases have been detected worldwide. From the total cases, only about 8 million cases ended in recovery, while the other 586,821 cases resulted in deaths. With a Case Fatality Rate (CFR) value of only 7 %, COVID-19 is not categorized as deadly as a MERS disease which has a CFR value of 34%.
There is an early report's advice that transmission might occur from three paths [2] i.e., (1) contact with an object that has the virus on it; (2) inhaling droplets emitted by sneeze and coughs; (3) micro-droplet infection. This third infection path is suspected to have resulted from a close contact (conversation or even from standing at a certain distance apart) with infected individuals.
In the experiment from NHK, it was found that when an infected individual coughs in a closed space, it is estimated that about a hundred thousand droplets are released within a few seconds.
Compared to large droplets that fall to the floor within 20-30 seconds, micro-droplets remain in the air for a longer period. This leaves people at risk of infection virtually all the time.
Several interventions have been put in place by the government to control the spread of COVID-19 such as (1) endemic prevention strategy (mitigation strategy) and (2) endemic reduction strategy (suppression strategy). This first strategy focused on slowing the spread of COVID-19 such as combining home isolation of suspect cases, social distancing, and so on. The second strategy aims to reverse the growth of COVID-19 through the improvement of medical treatment policy, research on medicine and vaccine for COVID-19, and so on [3]. These strategies 3 COST-EFFECTIVENSS AND BACKWARD BIFURCATION ON COVID- 19 MODEL have their own challenges whether it is from human behavior, minimum funding for prevention/control strategy, and so on.
Since COVID-19 was declared by WHO as a global pandemic, there have been many mathematical models introduced by authors to understand how COVID-19 spreads with the help of many aspects such as a model to understand the impact of non-pharmaceutical interventions (NPIs) [3], a model for understanding the effect of undetected cases [4,5], and a model for parameter estimation [6,7,8]. In this article, the result objective is focused on understanding the effect of micro-droplet transmission and its contribution to COVID-19 spread. Many reports have stated that medical resources are very limited in many countries. Therefore, we have included a saturated treatment rate in our model. The model is constructed with a basic SIS model combined with a new compartment for the droplets. Mathematical analysis regarding the stability analysis of the equilibrium points and the existence of backward bifurcation was conducted rigorously. For practical purposes for the COVID-19 control program, we constructed an optimal control problem based on our proposed model by adding three types of intervention as a time-dependent parameter: medical mask use, disinfectant intervention, and medical treatment. An optimal control problem analysis was conducted to determine the best strategy to reduce the number of COVID-19 infected individuals with an optimal intervention, which will depend on time. Some numerical experiments were conducted to provide a visual interpretation of the analytical results.
The rest of the paper is organized as follows. We propose our COVID-19 transmission models in Section 2. The first model is proposed without any intervention implemented. The purpose is to analyze the possible intervention strategy that might be implemented for COVID-19 prevention and endemic reduction purposes. The second model is proposed for policy purposes to eradicate COVID-19. In Section 3, optimal control simulations are conducted based on the second model.
Discussion on the cost-effectiveness analysis provided in Section 4, and followed with some conclusions in Section 5. To study the effect of indirect transmission and limited medical resources on the spread of COVID-19, in this section, we present the formulation of the model that will examine the dynamics of COVID-19. Human population ( ( )) is divided into two compartments: susceptible ( ( )) and infected ( ( )) compartments. The model is based on the standard SIS model. As mentioned previously, the droplet produced by an infected human through cough or sneeze can last a long time on surfaces [9]. Therefore, it is important to accommodate the existence of the droplets that are free in the environment which, in this article, we denote with ( ( )).

Model
Let the recruitment with the constant rate given by always be susceptible. The susceptible population is decreased by infection with an infected population following a population dependent infection function. This infection parameter is formed by a multiplication between the average number of meetings between and for each time interval Δ denoted as ℎ , and the chance of success infection denoted as ℎ . Therefore, the infection caused by direct contact between susceptible and infected individual is given by ℎ ℎ + , where + is the proportion of the susceptible population in time . For simplification, we wrote ℎ ℎ = ℎ . In addition to direct infection, we assumed that COVID-19 can also be transmitted to susceptible humans through contact between healthy human and surfaces containing the coronavirus . Unlike the previous direct infection, the indirect infection in this article is assumed to be modeled as a mass contact between the healthy population and coronavirus on the surface with an infection rate of .
Therefore, indirect transmission of Covid-19 is given by .
Next, it is important to mention other details on the recovery rate in this article. In many countries, medical resources significantly impact the success of COVID-19 control policy.
Unfortunately, the medical resources to control COVID-19 are limited in terms of the number of medical personnel, the availability of hospital beds, and many other factors. To accommodate this situation, the recovery rate that transfers the infected population into the susceptible population is in the term of 1+ , where and are the recovery rate and half saturated parameter, respectively.
Based on this explanation, we formulate the transmission model for COVID-19 with the direct and indirect transmission as well as saturated treatment as follows: where is the death rate induced by disease, is the average of production of coronavirus from the infected individual, and is the expected rate of survival of the coronavirus in the environment. System (1) is supplemented with a non-negative initial condition The rate of change on the total of human populations is as follows: As we assume that = ℎ , therefore we have the total of human population decreasing with respect to COVID-19 death incidence. It is easy to see that ( ) for all > 0, and ( ) ∈ [0, ℎ ] for → ∞. Since ( = 0) = , we have that ( ) always non-negative and has the upper bound ℎ . Therefore, we have ( ) ∈ [0, ℎ ] for all time > 0. With these properties, we have that system (1) is a biologically well-posed problem.

Model analysis
Before we analyze the model, it is important to reduce the model into a two-dimensional system, using the fact that has a very fast dynamic as the survival rate is much larger than the human survival rate ℎ , that is ≫ ℎ . Using the Quasi-Steady-State Approximation (QSSA) approach [10] and taking in the equilibrium state, we have * = . Substituting * into and on system (1), we have the new COVID-19 model: (2) The next step is re-scaling the time-space into recovery rate time-space by using = / , where is the new time-space. Thus, equation (2) is now read as follows: where * = , ℎ * = ℎ , * = , ℎ * = ℎ , and * = . Using these assumptions, all parameters now have no dimension except which is 1 . For the sake of written simplification, we neglect the "*" sign in all parameters and read system(3) as follows: One of the most important concepts in analyzing the epidemiological model of disease is the concept of the basic reproduction number which is commonly denoted as 0 . The basic reproduction number has been used by many authors to determine the behavior of the stability of the equilibrium from their disease model [11][12][13][14][15]. Commonly, the disease will die out whenever 0 < 1 and persist whenever 0 > 1. In this section, we will derive the basic reproduction number as the spectral radius of the next-generation matrix of the disease model. Please see [16] for further detail on the next-generation matrix approach in determining 0 and [17][18][19][20] for more examples of the implementation of this method in some epidemiological models. System (4) always has a trivial equilibrium state, that is, a COVID-19 free state which given by Γ 0 = ( , ) = ( ℎ , 0).
Applying the next-generation matrix approach to system (4), the basic reproduction number of system (4) is given by

COST-EFFECTIVENSS AND BACKWARD BIFURCATION ON COVID-19 MODEL
The eigenvalues of 0 are 1 = − ℎ and 2 = (1 + ℎ + )( 0 − 1). Therefore, it is clear that all eigenvalues of 0 will be negative if 0 < 1. Hence, Γ 0 is locally asymptotically stable if 0 < 1. On the other hand, if 0 > 1, then Γ 0 will be a saddle point. Hence, the proof is complete. ∎ From the form of the basic reproduction number, when = 0, which indicate that we neglect the indirect transmission, we have

Remark 1. The indirect transmission that occurs from contact with objects containing viruses on
their surface will increase the standard basic reproduction number from direct contact transmission only. Therefore, using disinfectant to kill the coronavirus could be a wise option to reduce the transmission of COVID-19 in the community.
Next, we analyze the existence of the untrivial equilibrium point. Taking the right-hand side of system (4) equal to zero and solving it with respect to and , the endemic state of system (4) is given by and 0 = ℎ ( + ℎ + 1)( 0 − 1), while 1 is taken from the positive roots of the following polynomial: while 2 has an expression too long to be shown in this article.
To analyze the existence of a positive solution of polynomial ( ), we will analyze the sign Substituting ℎ = ℎ * into ( ), we have * ( ) = 3 * 3 + 2 * 2 + 1 * + 0 * = 0, where 3 * , 2 * , 1 * and 0 * now is a function of 0 . Taking the derivative of respect to 0 from * ( ), and then solve it respect to It can be seen that 0 at = 0, 0 = 1 will be negative if Therefore, we can conclude that whenever > † , there will always exist an endemic state of COVID-19 of system (4) when 0 < 1. This result is summarized in the following theorem.

Theorem 2. System (4) will have a positive endemic state when
Theorem 2 states the existence of the COVID-19 endemic state even though 0 < 1. In the following theorem, we state the stability criteria of the COVID-19 equilibrium state.

Theorem 3. System (4) undergoes a backward bifurcation at
Proof. To analyze the stability of the COVID-19 endemic state, we will use the Castillo-Song theorem which introduced in [27]. First, let system (4) be redefined as follows Next, let ℎ be the bifurcation parameter. Hence, ℎ is evaluated when 0 = 1.
9 COST-EFFECTIVENSS AND BACKWARD BIFURCATION ON COVID-19 MODEL Substituting ℎ * , Γ 0 = ( ℎ , 0) into the Jacobian matrix of our system and evaluating the eigenvalues with respect to it, we have one zero simple eigenvalue, and the other eigenvalue is − ℎ . Hence, we can use the center-manifold theorem to analyze the stability of Γ 1 .

Hence, system (4) undergoes a backward bifurcation at
Hence, the proof is complete. ∎ Recall that is the saturation parameter for the treatment. Therefore, we have the following remark.
Remark 2. A more significant saturation parameter (indicated by lack of facility of treatment, the minimum number of beds in the hospital, minimum number of doctors, and other treatment facilities) will increase the probability of the existence of endemic equilibrium state of COVID-19, even though the disease-free state is stable when 0 < 1.

Numerical experiments
In this section, we illustrate the numerical examples of the possible phenomena of the COVID-19 model in (4), especially describing Theorem 2 numerically.
Using this parameter, we found that † = 10 −3 . Therefore, to conduct the forward bifurcation, we choose = 10 −4 < † . A bifurcation diagram for system (4) is shown in Fig.1 with respect to ℎ as the bifurcation parameter. In this setting, we found that system (4) only has the COVID-19 free equilibrium as the stable equilibrium whenever ℎ * < 1.0023 ( 0( ℎ =1.0023) = 1). When ℎ * increases, leaving ℎ * = 1.0023, the COVID-19-free state becomes unstable. On the other hand, the stable COVID-19 endemic starts to appear and increases when ℎ * increases. Biologically, this means that system (4) is robust with the initial condition as the long-term behavior of the model only depends on 0 whether its larger or smaller than 1. Example 2. Backward bifurcation. Similar to the previous example, from Section 3 we found that system (4) will have a backward bifurcation whenever > † . To illustrate this phenomenon, we use the same parameter value as example 1. Therefore, to conduct the backward bifurcation, we choose = 10 −2 > † . A bifurcation diagram for system (4) shows in Fig.2 with respect to ℎ * as the bifurcation parameter. In this setting, we found that system (4) only has the COVID-19 free equilibrium as the stable equilibrium whenever ℎ * < 0.385. When ℎ * = 0.385, hysteresis appears, which is shown by the appearance of another positive equilibrium.

Preliminary conclusions
In this section, we propose an epidemic model with a saturation recovery rate from infected to susceptible population considering the indirect transmission from the free pathogen. Our mathematical analysis started by approximating the coronavirus compartment in its equilibrium condition using the QSSA approach. Our model has two equilibrium states. The COVID-19 free state is locally stable if the basic reproduction number is less than one and becomes unstable otherwise. Our analysis shows that it is possible that the COVID-19 endemic state might appear even if the basic reproduction number is less than unity. This condition might lead to a misinterpretation of the condition of the COVID-19 eradication program by the policymaker in the community. Our analysis also shows that the indirect transmission which occurs from touching objects that contain viruses on their surface could increase the possibility of the existence of COVID-19 in the community. Therefore, except only using social/physical distancing and rapid assessment [21] or medical treatment [22], the use of disinfectant to wipe out the coronavirus that could potentially cause indirect transmission could be considered an additional effort to reduce COVID-19 transmission. Hence, we will improve our proposed model to consider three different interventions: medical mask use, disinfectant, and treatment. These three interventions will be introduced as a time-dependent parameter in the following section.

Improved model with control program
Based on our previous analysis, it is important to consider other interventions except focusing on social distancing to eradicate COVID-19. Hence, in this section, we modify our previous COVID-19 transmission model in (1). The modification entails adding three different interventions with specific purposes; they are as follows: (i) Medical mask. The medical mask is strongly recommended by many policymakers in many countries not only for the infected individual but also for the susceptible one as many infected people do not show any symptoms. In our model, this intervention is denoted by 1 ( ). The construction of the model is as follows. First, we assume that only 1 ( ) proportion of the human population use a medical mask, while the rest (1 − humans from direct contact with infected individuals perfectly with the efficacy of the medical mask being 100%. Therefore, the infection caused by contact with an infected individual is 0. Hence, infection from direct contact only occurs from contact between a susceptible and an infected individual who does not wear any medical mask, i.e. (1-1 ( )) ( ) and (1-1 ( )) ( ), respectively. Therefore, the new infection coming from direct contact with an infected individual is given as follows: (ii) Disinfectants. Disinfectants are implemented to kill coronavirus, which probably lies in the surface that already touch/had droplet by an infected individual previously. In our model, we denote this intervention as 2 ( ).
(iii) Treatment. Although specific medicine to cure infected individuals has not yet been created, medical support can still increase the recovery rate of an infected individual.
We use 3 ( ) to symbolize this intervention on our model.

The basic reproduction number
Assuming all interventions are constant in time, system (6) has a COVID-19 free equilibrium given by Using the next-generation matrix method [16], the basic reproduction number of system (6) is given by 14 DIPO ALDILA 0 = + = ℎ (1 − 1 ) 2 ( ℎ + + + 3 ) + ℎ ( + 2 )( ℎ + + + 3 ) . (8) Please note that denotes the basic reproduction number for the direct transmission caused by contact between the susceptible and the infected individual, while is the infection from the indirect infection. From [24] we have the following theorem.

Characterization of the optimal control problem
In this section, we analyze the optimal control problem for COVID-19 transmission model (6) to determine the optimal trajectories of ( ), ( ), and ( ) which correspond to the optimal intervention strategy. Let us define the objective function as follows: ( 1 , 2 , 3 ) = ∫ ( 1 ( ) + 2 ( ) + 1 2 1 1 2 + 1 2 2 2 2 + 1 2 3 3 2 ) . 0 (9) This objective function describes the aim of our optimal control problem: minimizing the number of infected people and free coronavirus in the environment with a minimum cost of intervention strategies. We assume that the cost of the treatment is quadratic as has already been used by many authors in [12,17,20,24]. The positive coefficients 1 , 2 , 1 , 2 , and 3 are the weight parameters which will guarantee the balance of each component in the cost function (9).
The Pontryagin Maximum Principle [25] has been used in this article to derive the necessary conditions that the optimal control parameters should satisfy. According to this principle, 1 * , 2 * , and 3 * as the optimal trajectories of control parameters with corresponding optimal * , * , and * will minimize the cost function ( 1 , 2 , 3 ) for a fixed final time .
Proof. The co-state system in (16) is the direct result from equation (12) in which we differentiate the Hamiltonian with respect to each state variable. The transversality condition is derived from (15). Taking the derivative of with respect to 1 , solving it with respect to zero, and combining it with its lower bound and upper bound 0 and 1 respectively, we have 1 * ( ) in equation (17). A similar approach is used to find 2 * and 3 * . ∎

RESULTS IN THE OPTIMAL CONTROL PROBLEM
To conduct the numerical simulation for the optimal control problem in this article, we use the socalled forward-backward iterative method [26] and MATLAB codes for the implementation.
Generally, the forward-backward iterative method is as follows. We begin with having an initial estimation for all control variables, solving the state variables in system (6) forward in time for all time ∈ [0, ], and then calculating the value of initial cost function (9). Having the solutions of all state variables, we use the same to find the solution to all co-state variables using system (16) backward in time. We update our control variables using (17). Using these updated control variables, we re-calculate the state variables and cost function. We repeat this step until some convergence condition is achieved, which in this case is that the norm-2 difference of control and state variables in iteration + 1 compared to the ℎ iteration should be less than small number .
To conduct the simulation in this section, we use the following parameter values: In simulating the optimal control problem, we use three different strategies, i.e. the use of single intervention (Section 3.1), two interventions (Section 3.2), and all interventions (Section 17 COST-EFFECTIVENSS AND BACKWARD BIFURCATION ON COVID-19 MODEL 3.3). To run the optimal control simulation, we use the aforementioned parameter value, weight parameter as follows 1 = 1, 2 = 1, 1 = 10 12 , 2 = 10 12 , 3 = 10 17 , to balance the cost function and describe the cost for implementation. It is evident that 1 = 2 as the cost for medical mask use is almost the same as the cost for personal disinfectant use. On the other hand, 3 is higher than 1 and 2 as the cost for medical intervention requires a more expensive cost for implementation. Furthermore, we use the following state initial condition: ( 0 , 0 , 0 ) = (99.000,1.000,0).

Single intervention
It can be seen from Fig. 3-5 that without the implementation of any control, the dynamic of infected peoples increased to outbreak at Day 40 and started decreasing to equilibrium when started increasing. Consequently, the implementation of control will behave related to the dynamic of the infected individual. It can be seen that 1 and 2 are maximum from the beginning of the simulation time until a few days after the outbreak; they then start to decrease to 0. On the other hand, the dynamic of 3 is monotonically decreasing from the beginning of simulation time. This is probably because the cost of treatment in the hospital is high compared to other interventions.

Combination of two interventions
In Fig. 6, we can see that each strategy could reduce the number of predicted infected individuals with a combination of 1 and 2 to yield a better result (please see the green curve in Fig. 6) with the cost function = 6.844 × 10 11 . Second-best interventions are provided by a combination of 1 and 3 (please see the blue curve in Fig. 7) with the cost function = 3.945 × 10 14 . Compared to 1 and 2 only, the number of infected individuals increases in the time close to the final time of simulation which indicates that possible delayed outbreak will occur.
From 8, a combination of 2 and 3 should be given at a maximum rate at the beginning of the simulation and start to decrease when the peak of the outbreak passes.

Combination of all interventions
From Fig.9, it is evident that the number of infected individuals can be reduced significantly from the beginning of simulation time to be less than 1000 individuals for the entire simulation. For the control interventions, it can also be seen that they are only present at a high rate before the predicted outbreak (without control strategy) occurs; all controls start to decrease to 0 after that. It is also apparent that 1 maintains the highest rate longer than other controls. The cost for intervention when all interventions are in place is 4.9 × 10 13 .

COST-EFFECTIVENESS ANALYSIS
In this section, we perform a cost-effectiveness analysis of our optimal control Strategy discussed in Section 3. We consider two types of methods to determine the best strategies, namely the infection averted ratio (IAR) approach and the average cost-effectiveness ratio (ACER) approach.
The IAR is defined as = number of infection averted Number of recovered .
The highest IAR ratio is the best strategy. On the other hand, = Total cost produced by intervention Total number of infection averted .
The smallest value of ACER is the best strategy. The result for numerical calculation according to total cost function (9), number of infected averted, total recovery, cost for intervention, IAR and ACER can be seen in Table 1

Infected averted ratio (IAR)
According to Fig. 9, the best strategy is Strategy 7 when all controls are implemented. The second and third best are Strategy 4 (when medical mask and disinfectant combined) and Strategy 5 (when the medical mask and medical treatment are combined), respectively. We can see that Strategy 3 is the least cost-effective.

DISCUSSION AND CONCLUDING REMARKS
In this study, a deterministic model for understanding the spread of COVID-19 that considers direct and indirect transmission has been proposed as a system of ordinary differential equations.
Mathematical analysis of the existence and local stability criteria of all equilibrium points was conducted and indicated how the basic reproduction number plays an important role in determining the long-term behavior of the model. Our proposed model might exhibit a backward bifurcation phenomenon at 0 = 1 if the saturated treatment parameter is larger than its threshold.
Further improvement of the model under the interventions medical mask, disinfectant, and medical treatment were proposed later. The model has a disease-free equilibrium which is locally stable when the basic reproduction number is less than one. The optimal control problem of the proposed model constructed is based on Pontryagin's Maximum Principle to minimize the number of infected individuals with an optimal cost. The problem was solved using the forward-backward iterative method.
The numerical simulation of the optimal control problem was conducted using three strategies.
In the first Strategy, when only a single intervention was implemented, our simulation revealed that the medical mask use is the best single intervention strategy. This was followed by disinfectant and medical treatment. In the second Strategy, a combination of medical mask use and disinfectant yielded a better result compared to the other combination. In the last Strategy, a combination of all 26 DIPO ALDILA interventions together provided the best result in reducing the number of infected individuals compared to the best strategy in the first and second strategies. All strategies show the profile of control adapted with the number of infected individuals in the community. A larger intervention is needed to substantially reduce the number of infected individuals, and this intervention then start to reduce whenever the number of infected individual start to show a trend to decrease. Further analysis to determine the most cost-effective strategy was conducted using the IAR and ACER indices. We found that using all controls simultaneously is the most cost-effective method based on the calculation of IAR and that using a medical mask only is the most cost-effective method in a single intervention which was analyzed using the ACER calculations.
The results in this study provide an alternative scientific background for the eradication of COVID-19 in the community. The limitation of budget that will force policymakers to save on intervention costs by choosing priority strategies has been discussed here. Each intervention has its own focus. If the policy is to reduce the number of deaths, then the intervention should focus on medical treatment. However, this intervention comes with high-cost consequences. On the other hand, if the policy is to prevent the outbreak of COVID-19, then medical mask use or disinfectant intervention or combination of both is the best strategy as that could prevent the outbreak well. It also costs less compared to medical treatment intervention. For future work, including other interventions into the model could be considered such as physical/social distancing intervention, rapid infected assessment, quarantine, and many more.