Adv Syst Sci Appl 2022; 01:148–166 Published online at https://ijassa.ipu.ru. Sensitivity and Optimal Control Analysis of Japanese Encephalitis Disease: A Mathematical Model Naba Kumar Goswami1* 1Department of Mathematics, PET Research Foundation, PES College of Engineering, Karnataka, India Abstract: Japanese encephalitis is caused by flavivirus, which can affect both humans and animals. A Japanese encephalitis disease model is proposed and analyzed to investigate the transmission dynamics of Japanese encephalitis with saturated incident rate and saturated treatment function. Assume that the infected pig populations are growing logistically in the environment. The basic reproduction number (R0) of the model is obtained using the next- generation matrix method. The sensitivity analysis is performed to identify the key parameters that affect the basic reproduction number, which can be regulated to control the transmission dynamics of the disease. The model is extended to optimal control model incorporating three time- dependent inputs for the control of transmission route by using Pontryagin’s maximum principle. This study is significant as Japanese encephalitis disease poses serious challenges for public health in past few years. Numerical simulation is performed to support our analytical findings. Keywords: japanese encephalitis, basic reproduction number, sensitivity, optimal control, simulation 1. INTRODUCTION The Japanese encephalitis is primarily a mosquito-borne rural disease, which is transmitted to human through the bite of an infected Culex species mosquitoes. Japanese encephalitis virus (JEV) is a flavivitus of the family of flaviviridae like the viruses of the Yellow fever, West Nile, and Dengue. The virus is transmitted in an enzootic cycle among mosquitoes and vertebrate amplifying hosts, primarily domestic pigs and/or water birds (Wild Heron) [1, 2]. Most impotently there is no such type of evidence for human to human transmission. The first case of JEV disease was recorded in 1871 in Japan [3]. The clinical cases of Japanese encephalitis globally estimated every year nearly 68,000, with approximately 13,600 to 20,400 deaths. The Japanese encephalitis primarily affects children less than 15 years of age. The fatality rate varies between 20%− 40%, but it may reach 58% and over. The incubation period for Japanese encephalitis is between 4–14 days [1, 3]. The primary symptoms of the disease includes fever, chills, headache, and vomiting. Other symptoms like neurological symptoms, drowsiness, dizziness, confusion, abdominal pain and diarrhea might develop over the next few days [5]. There is no specific medicine for treatment, but patients should be hospitalized for supporting care. During the infected period of Japanese encephalitis, people should drink plenty of fluids and rest sufficiently. The vaccine is available to control Japanese encephalitis virus disease but it does not protect 100 %. Japanese encephalitis virus has been identified most common in rural areas in southeast Asia, Eastern Asia, Southeastern Asia, and Pascific Island. The disease is investigated globally in the following territories : Australia, Bangladesh, Bruni, Burma, Cambodia, China, Guam, India, Indonesia, ∗Corresponding author: nabakrgoswami@gmail.com SENSITIVITY AND OPTIMAL CONTROL ANALYSIS OF JAPANESE 149 Japan, Laos, Malaysia, Nepal, North Korea, Pakistan, Phillipines, Russia, Saipan, Singapore, South Korea, Sri Lanka, Taiwan, Thailand, Timor-Leste, Vietnam. The main vector of Japanese encephalitis is C. Annulirostris. In Indian condition many other secondary vectors have reported including Anopheles subpictus, A. peditaeniatys, C.epidesmus, C.gelidus, C.whitmorei, C. pseudovishnui, M. uniform and Mansonia indiana. Pigs are considered to be principal amplifying host in India. The transmission is possible between May and September. The first case of Japanese encephalitis disease was diagnosed in 1955 in Tamil Nadu in India [4]. Consequently, several major out-breaks have been listed from different parts of the country, mostly in rural areas of Assam, Bihir, Harayana,Karnataka, Tamil Nadu, Uttar Pradesh. In India, 1,500 to 4000 JEV cases are reported every year [3]. According to National health mission in Assam most of the upper Assam districts, including Dibrugarh, Golaghat, Jorhat, Lakhimpur and lower Assam district of Kamrup have been badly affected every year by the Japanese encephalitis mosquito-borne disease. Several mathematical models have been studied and analyzed the dynamics of Japanese encephalitis transmission [3-13] as yet. In [3], authors discussed the current knowledge of the epidemiology and the parthenogenesis of this deadly disease have been summarized. In [4], authors deals with application of various control strategies to Japanese encephalitis among human, pig and mosquito. In [5], authors proposed a mathematical model for the spread of Japanese encephalitis with emphasis on the environmental effects on the aquatic phase of mosquitoes. In [6], authors proposed a mathematical model with impact of vaccination on of JE with standard incidence rate of mosquitoes, pigs and humans. In [7], authors has been studied stability and bifurcation analysis of Japanese encephalitis model with/without effect of some control parameters. In [8] authors developed a mathematical model and analysis the spread of Japanese encephalitis with environmental effects. In [9] authors build an SIRS epidemic model of Japanese encephalitis and discussed analytically. In [10], authors have formulated a deterministic model with saturated treatment function and found that north east states of India need better treatment and awareness about the disease. In [11], author contracted a mathematical model with saturated incidence function and analysis optimal control theory. In [12], author discussed different causes of backward bifurcations in some epidemiological models. In [13] authors formulated a mathematical model and regulate important parameters in the spread of the disease through the sensitivity analysis. In this paper, we have constructed a Japanese encephalitis mathematical model by considering standard incidence type interaction with saturated incident and saturated treatment. As India is a developing countries, so treatment may not available for all if the number of infectives are very large. So in the present study we have incorporated the saturated treatment as discussed in [10]. The remaining part of the paper is organized as follows: Section 2 describes the model and established biological feasible region; Section 3 presents basic reproduction number; Section 4 deals with data; Section 5 deals with sensitivity analysis of basic reproduction number R0; Section 6 describes optimal control problem; section 7 demonstrates the numerical simulation of optimal control and finally Section 8 concludes the paper. 2. THE MODEL FORMULATION The model consists six different compartments such as Susceptible humans (Sh), Infected humans (Ih) and Recovered humans (Rh), Susceptible Japanese encephalitis mosquitoes (Sj), Infected Japanese encephalitis mosquitoes (Ij), and Infected pigs (Ip). Many epidemic models has been used bilinear incident rate to prevent and control the spread of the infectious disease. Here we introduced saturated incident rate in the proposed model as saturated incident rate is more suitable than bilinear incident rate. Based on following assumption the system of differential equations formulated. 150 N.K. GOSWAMI • β1ShIj 1 + kIh is the saturated incident rate, which tends to a saturated level when Ih gets large, β1Ij measures the infection force when disease is entering a fully susceptible population [14]. • β1ShIj is known as bilinear incidence rate and 1 1 + kIh measures the inhibition effect from the behavioral change of susceptible individuals when their number increases or from the crowding effect of the infective individuals [9, 14, 15]. • 1 1 + kIh is the incidence rate which is more reasonable than the β1ShIj (bilinear incidence rate), [10]. • Assume that through treatment the infected individuals recover at a saturated function h(Ih) = γhIh 1 + αIh , where γh is recovery rate. When α = 0, the saturated treatment function becomes linear [14, 16–18]. • The model is a non-linear standard incidence type interaction between human and mosquitoes. • The individuates are recruited in the region at a constant rate Λh and they join the susceptible class. • The susceptible individuals become infected when the mosquito bites and join the infected human class at the rate of β1. • Infected individuals get treatment and recover from encephalitis with recovery rate γh and join recover class. • The susceptible Japanese encephalitis mosquitoes becomes infected and join infected mosquitoes class at the rate β2. • Susceptible mosquitoes bite the infected pigs and becomes infected at the rate β3. • The infected pig population are growing logistically. Keeping in view the above facts, the mathematical model is framed as follows: dSh dt = Λh − β1Sh ( Ij 1 + kIh ) − µhSh + δhRh dIh dt = β1Sh ( Ij 1 + kIh ) − γhIh 1 + αIh − (µh + µ1)Ih dRh dt = γhIh 1 + αIh − (δh + µh)Rh dSj dt = Λj − β2Sj Ih Nh − β3SjIp − µjSj dIj dt = β2Sj Ih Nh + β3SjIp − µjIj dIp dt = rIp ( 1− Ip K ) − µpIp (2.1) 2.1. Positivity and Boundedness Here, we shall show that the positivity and boundedness of the population. From the system (1), we have dSh dt ∣∣∣ Sh=0 = Λh + δhRh > 0, dIh dt ∣∣∣ Ih=0 = β1ShIj ≥ 0, SENSITIVITY AND OPTIMAL CONTROL ANALYSIS OF JAPANESE 151 Table 2.1. Description of biological parameters Parameter Description Value Λh Rate of recruitment of human population 20 Λj Rate of recruitment of vector(mosquito) 40 β1 Transmission probability interaction between Sh and Ij 0.00256 β2 Transmission probability interaction between Sj and Ih .031 β3 Transmission probability interaction between Sj and Ip 0.013 α A positive constant 0.02 γh Recovery rate of human population 0.005 δh Rate at which recovered human become susceptible due to loos of immunity 0.05 µh Natural mortality rate of human population 0.0421 µ1 Mortality rate of human due to infection 0.0326 µj Natural mortality rate of vector 0.112 µp Natural mortality rate of pigs 0.51 r The growth rate coefficient of infected pigs 2 K Carrying capacity of mosquito population 10 dRh dt ∣∣∣ Rh=0 = γhIh 1 + αIh ≥ 0, dSj dt ∣∣∣ Sj=0 = Λj > 0, dIj dt ∣∣∣ Ij=0 = β2Sj Ih Nh + β3SjIp ≥ 0, dIp dt ∣∣∣ Ip=0 ≥ 0 Here, all the rates are non-negative, so if we start in the interior of the non-negative bounding R6, we shall always remain in this cone keeping mind of the fact that direction of the vector field is inward on all the bounding planes. We note the change rate of the total population Nh = Sh + Ih +Rh and Nj = Sj + Ij are given by the following differential equations : dNh dt = Λh − µhNh − µ1Ih, dNj dt = Λj − µjNj. This gives lim sup t→∞ Nh ≤ Λh µh , lim sup t→∞ Nj ≤ Λj µj . Therefore, all Sh(t), Ih(t), Rh(t) are bounded by Λh µh and the solutions Sj(t), Ij(t) are bounded by Λj µj . Hence, the biological feasible region of the proposed system (1) is given by the following positively invariant region: Ω = (Sh, Ih, Rh, Sj, Ij, Ip) ∈ R6 + : (Sh + Ih +Rh) ≤ Λh µh , (Sj + Ij) ≤ Λj µj . 152 N.K. GOSWAMI 3. BASIC REPRODUCTION NUMBER The disease-free equilibrium for the system (2) as follows E0=(N0 h , I 0 h, R 0 h, N 0 j , I 0 j , I 0 p ) =( Λh µh , 0, 0, Λj µj , 0, 0 ) . The reproduction number R0 gives the average number of infected individuals generated by the one in a fully susceptible population and for our model, it is given by the above expression of R0. Mathematically, R0 is a threshold parameter for the stability of a disease- free equilibrium and is related to the highest and final size of an epidemic. To find the basic reproduction number R0, we consider ( dIh dt , dRh dt , dIj dt , dIp dt ) and using the next-generation matrix method as described in [19, 20]. The matrix F and V as follows: F =  β1(Nh − Ih −Rh) Ij 1 + kIh 0 β2(Nj − Ij) Ih Nh + β3(Nj − Ij)Ip 0  , V =  γhIh 1 + αIh + (µh + µ1)Ih − γhIh 1 + αIh + (δh + µh)Rh µjIj −rIp ( 1− Ip K ) + µpIp  Now the matrix F ( Jacobian of F at disease-free equlibrium E0) and V (Jacobian of V at disease-free equilibrium E0) are obtain as F =  0 0 β1N 0 h 0 0 0 0 0 β2 N0 j N0 h 0 0 β3N 0 j 0 0 0 0  and V =  γh + µh + µ1 0 0 0 −γh µh + µp 0 0 0 0 µj 0 0 0 0 −r + µp  The largest eigenvalue of FV −1 is called the basic reproduction number R0 and is obtained as follows: FV −1 =  0 0 β1N 0 h µj 0 0 0 0 0 β2 N0 j N0 h(γh + µh + µ1) 0 0 β3N 0 j µp − r 0 0 0 0  , SENSITIVITY AND OPTIMAL CONTROL ANALYSIS OF JAPANESE 153 2005 2010 2015 2020 0 200 400 600 800 Po pu la tio ns Years Cases Deaths Fig. 4.1. Trend of Japanese Encephalitis cases in Assam, India during 2004-2019 The spectral radius of the matrix FV −1 is called the basic reproduction number R0 and is obtained as follows: R0 = √ β1β2Λj µ2 j(γh + µh + µ1) The quantity R0 is knownas basic reproduction number, the expected number of secondary cases produced in completely susceptible population, by a typical infective individual, 4. DATA SCENARIO In India, approximately 597,542,000 people are live in Japanese encephalitis region, and every year 1,500 to 4,000 japanese encephalitis cases are reported. Here, from 2004 to 2019 annual reported cases of japanese encephalitis are using for graphical representation. For this study, we consider two states of India which are Assam and Uttar Pradesh as these are Japanese Encephalitis infected states of India. The total reported cases from 2004 to 2019 in Assam and Uttar Pradesh are 7,325 and 33,572, which are demonstrated in Figure 1 and 2 respectively. In Uttar Pradesh highest number of cases were reported in 2006 and lowest number of cases were reported in 2019. In Assam highest number of cases were reported in 2014 and lowest number of cases were reported in 2005. 154 N.K. GOSWAMI 2005 2010 2015 2020 0 2000 4000 6000 Po pu la tio ns Years Cases Death Fig. 4.2. Trend of Japanese Encephalitis cases in Uttar Pradesh, India during 2004-2019 5. SENSITIVITY ANALYSIS OFR0 To visualize the effect of the parameters involved in the expression of basic reproduction number R0, we calculated the normalized forward sensitivity indices of R0 to those parameters [21]. We keep the parameters at the same values as in Table 1. The normalized forward sensitivity index of a variable to a parameter is determined by the ratio of the relative change in the variable to the relative change in the parameter [13]. The sensitivity indices of R0 with respect to the parameters of interest are depicted in Figure 3, 4 and 5. The figure shows that when the parameters β1, β2, and Λj increase, keeping the other parameters constant, the value ofR0 increases as these parameters have positive indices. Instead, increase in the values of parameters γh, µh, µ1, and µj leads to decrease in the values of R0 as they have negative indices. It is noted that the sensitivity index of R0 is 0.5 for the parameters β1, β2,Λj . It means that 0.5% increase in the values of any of these parameters, keeping other parameters fixed, will result in 0.5% increment in the value of R0. Precisely, lower values of R0 are of crucial importance as they increase the chances of disease eradication. Therefore,it is imperative to prevent an increase in the parameters having positive sensitivity indices whereas increasing the values of parameters having negative indices is instead preferred.Thus, any prevention measure aiming at reducing the former parameters and increasing the latter must be taken into serious consideration. If such variable is differentiable with respect to the parameter, then the sensitivity index is defined using partial derivatives, [13]. The normalized forward sensitivity index of R0, which is differentiable with respect to a given parameter P, is defined by Y R0 P = P R0 ∂R0 ∂P SENSITIVITY AND OPTIMAL CONTROL ANALYSIS OF JAPANESE 155 1 2 h h 1 j j Parameters -1 -0.5 0 0.5 S en ci tiv ity In di ce s Fig. 5.3. Normalized forward sensitivity indices of R0 with respect to model parameters. Parameter values : β1 = 0.00256, β2 = 0.031, γh = 0.245, µh = 0.0421, µ1 = 0.0326, µj = 0.112,Λj = 40 The above formula can be used to compute the analytical expression for the sensitivity of R0 to each parameter that it includes. Plotted contour plot of some key parameters on R0. If we increase or decrease the value of β1 and β2, and then the value of R0 will increase or decrease, respectively, which is shown in Figure 4(a). Also, we can conclude that increasing the value of Λj will increase the value of the R0, but if we increase the value of µj , it will decrease the value of the R0, which is the death rate of the mosquito, and this is the one of the best control policy to reduce the disease infection from the population, manifest in Figure 4(b). Similarly other significant changes ofR0 are shown in Figures 4(c) and 5(a,b,c) based respective parameters. So, the correct estimation of these parameters is very important to predict transmission of this disease. 6. OPTIMAL CONTROL PROBLEM In this section, three different types of control intervention viz., u1(t), u2(t), u3(t) are incorporated into the model system (1) and extended to optimal control problem. There interventions are implementing either pharmaceutical(treatment) or non- pharmaceutical(effect of information). The main goal this research is to investigate the best control strategies with minimum cost of implementation as well as financial loss generated. However, the vaccination and effective medicines are the most useful strategies to prevent contagious transmission. The details of each intervention are described as follows: • Control variable u1(t): The force of Japanese Encephalitis infections is reduce by (1- u1(t)), where u1(t) measures effect due to the use of electronic devices, insecticide- treated bed nets, and mosquito repulsive lotions are used to reduce transmission between mosquitoes and human. If one can reduce between mosquitoes to human transmission rate, then it will be helpful in controlling the mosquito-borne disease. Hence, we consider u1(t) as a control intervention to reduce human-mosquito interaction. • Control variable u2(t): The control variable u2(t) is the pig vaccination lowered the mosquito infection rate. Pig vaccination is associated with a reduction in human cases. In this an effective intervention but requires considerable amount of money. Thus, our main objective is to minimize the disease effect in pig population with minimum cost. • Control variable u3(t): The control variable u3(t) is the use of effective medicines for the treatment measure of infectious humans, proper vaccination and insecticide to the Japanese Encephalitis infection causing species [22]. η is the modification parameter 156 N.K. GOSWAMI 0 0.2 0.4 0.6 0.8 1 1 0 0.5 1 2 10 20 30 40 (a) 0 0.2 0.4 0.6 0.8 1 j 0 0.2 0.4 0.6 0.8 1 j 0.01 0.02 0.03 0.04 (b) 0 0.2 0.4 0.6 0.8 1 h 0 0.2 0.4 0.6 0.8 1 1 6 7 8 9 10 (c) Fig. 5.4. Contour plots of the basic reproduction numberR0 with respect to (a) β1 and β2 (b) Λj and µj , (c) γh and µ1 SENSITIVITY AND OPTIMAL CONTROL ANALYSIS OF JAPANESE 157 0 0.2 0.4 0.6 0.8 1 1 0 0.2 0.4 0.6 0.8 1 h 5 10 15 20 25 (a) 0 0.2 0.4 0.6 0.8 1 1 0 0.5 1 h 5 10 15 (b) 0 0.2 0.4 0.6 0.8 1 2 0 0.2 0.4 0.6 0.8 1 j 0.1 0.2 0.3 0.4 0.5 0.6 0.7 (c) Fig. 5.5. Contour plots of the basic reproduction numberR0 with respect to (a) β1 and µh, (b) β1 and γh, (c) β2 and Λj 158 N.K. GOSWAMI in control u3(t). In this an effective intervention but requires considerable amount of money. Thus, our main objective is to minimize the disease effect in human population with minimum cost. Based on different health care restriction , it is essential to impose some bounds on controls as 0 ≤ u1(t), u2(t), u3(t) ≤ 1. If u1(t), u2(t), and u3(t) are equal to zero, then there is no effort being placed in these controls at time t and if they are equal to one then the maximum effort is applied. Keeping the view of the above assumptions, the optimal control model is formulated as follows: dSh dt = Λh − (1− u1(t))β1Sh ( Ij 1 + kIh ) − µhSh + δhRh dIh dt = (1− u1(t))β1Sh ( Ij 1 + kIh ) − ηu3(t)γhIh 1 + αIh − (µh + µ1)Ih dRh dt = ηu3(t)γhIh 1 + αIh − (δh + µh)Rh dSj dt = Λj − (1− u1(t))β2Sj Ih Nh − (1− u2(t))β3SjIp − µjSj dIj dt = (1− u1(t))β2Sj Ih Nh + (1− u2(t))β3SjIp − µjIj dIp dt = rIp ( 1− Ip K ) − µpIp (6.2) Consider u1(t) = u1, u2(t) = u2, and u3(t) = u3 for further analysis. 6.1. Cost Construction and Characterization of Optima Controls This particular section devotes into two parts which include determination of the total cost generated due to applied controls as well as diseases itself. Whereas the second part determines the analytical forms of the controls. 6.1.1. Total Cost Determination Here, the total cost is determined for the applied control interventions and disease as well, which need to be minimized. • Cost due to disease: The cumulative cost incurred due to the disease burden is modeled as follows: ∫ T 0 (X1Ih +X2(Sj + Ij) +X3Ip)dt This cost consists of various components such as cost due to loss of manpower, opportunity loss and other related wealth loss. • Cost incurred in insecticide-treated bed nets: The total cost involved in reducing the transmission between human to mosquito is given as∫ T 0 Y1u 2 1dt. • Cost incurred intreatment and vaccination for pig population: The total cost involved in reducing pig infection is defined as∫ T 0 Y2u 2 2dt. SENSITIVITY AND OPTIMAL CONTROL ANALYSIS OF JAPANESE 159 • Cost incurred in treatment for infected human population: The cumulative cost in process of treating infected human is defined by∫ T 0 Y3u 2 3dt. The cost functional corresponding to total cost incurred, for fixed time T , which need to be minimized is given by J(u1, u2, u3) = ∫ T 0 [ X1Ih +X2(Sj + Ij) +X3Ip + 1 2 (Y1u 2 1 + Y2u 2 2 + Y3u 2 3) ] dt (6.3) subject to the model system (2). The parameterX1 ≥ 0,X2 ≥ 0, Y1 ≥ 0, Y2 ≥ 0, and Y3 ≥ 0, are the weight and balancing constants, which measure the respective cost involvement over the interval [0, T ]. The term Y1u1 2 denotes the cost associated with Insecticide-treated bed nets are a form of personal protection that has been shown to reduce mosquito borne disease i.e. zika, dengue, malaria etc. and death due to mosquito borne disease in endemic regions. But at the same time, the execution of ITNs requires lot of money which is given by the term Y1u 1 1 [11]. Also, the total cost incurred in process primarily includes the cost involved in awareness and educational campaigns, cost of required manpower etc. The term Y2u2 2 signifies the cost associated with pig vaccination to reduce the mosquito infection rate and the term Y3u3 2 represents the cost associated with treatment for infected human population which primarily medical facilities, Hospitals expenditure, manpower, human vaccination etc [21]. In order to find an optimal control, u1∗,u2∗, and u3∗ such that J(u1 ∗, u2 ∗, u3 ∗) = min (u1,u2,u3)∈U J(u1, u2, u3), (6.4) where U is the control set and is defined as U = {(u1, u2, u3) : 0 ≤ u1, u2, u3 ≤ 1, t ∈ [0, T ]} Here, all the controls are bounded and measurable. 6.1.2. Existence and Characterization of Optima Controls Here, we shall first establish the existence of such control functions that minimises the cost functional J . The Lagrangian L of this problem is defined as: L(Ih, Sj, Ij, Ip, u1, u2, u3) = X1Ih +X2(Sj + Ij) +X3Ip + 1 2 Y1u1 2 + 1 2 Y2u 2 2 + 1 2 Y3u 2 3 Now, we shall use Pontryagin’s maximum principle [?] for necessary conditions for optimal controls system (2). For that by choosing A = (Sh, Ih, Rh, Sj, Ij, Ip), U = (u1, u2, u3) and λ = (λ1, λ2, λ3, λ4, λ5, λ6), the associated HamiltonianH can be written as H(A,U , λ) = L(Ih, Sj, Ij, Ip, u1, u2, u3) + λ1 dSh dt + λ2 dIh dt + λ3 dRh dt + λ4 dSj dt + λ5 dIj dt + λ6 dIp dt (6.5) Since u∗1,u ∗ 2, u∗3 are solutions to the control problem (2), there exists the adjoint variables λ1, λ2, λ3, λ4, λ5, λ6 satisfying the following conditions. dx dt = ∂H(t, x, u∗1, u ∗ 2, u ∗ 3, λ1, λ2, λ3.λ4, λ5, λ6) ∂λ 0 = ∂H(t, x, u∗1, u ∗ 2, u ∗ 3, λ1, λ2, λ3.λ4, λ5, λ6) ∂u dλ dt = −∂H(t, x, u∗1, u ∗ 2, u ∗ 3, λ1, λ2, λ3.λ4, λ5, λ6) ∂x (6.6) 160 N.K. GOSWAMI Theorem 6.1: For the objective functional (3) and the control set (6) subject to control system (2) there exists an optimal control u∗ = (u1 ∗, u2 ∗, u3 ∗) ∈ U such that J(u1 ∗, u2 ∗, u3 ∗) = min U J(u1, u2, u3). Proof To establish this result, we follow the Theorem 4.1 mentioned in [23] for the existence of optimal controls. As, we have discussed above that all the state variables (population) are bounded for each bounded controls coming from the control set U . Furthermore, Lipschitz condition with respect to state variables is satisfied by the right hand part of the model system (2). The control variable set U is also convex and closed by the definition and the model system (2) is linear in control variables. The integrand of the function L = X1Ih +X2(Sj + Ij) +X3Ip + 1 2 Y1u1 2 + 1 2 Y2u 2 2 + 1 2 Y3u 2 3 is convex on the control set U due to quadratic nature of control variables. Moreover, L = X1Ih +X2(Sj + Ij) +X3Ip + 1 2 Y1u1 2 + 1 2 Y2u 2 2 + 1 2 Y3u 2 3 ≥ 1 2 Y1u1 2 + 1 2 Y2u 2 2 + 1 2 Y3u 2 3. Now consider c1 = min(X1, X2, X3) > 0 and g(u1, u2, u3) = c1(u1 2 + u2 2 + u3 2). Thus, L ≥ g(u1, u2, u3) holds true and g is continuous. Also, g satisfies the condition |(u1, u2, u3)|−1g(u1, u2, u3)→ ∞ whenever |(u1, u2, u3)| → ∞. Thus, all the conditions for the existence of controls are fulfilled (for more details one can follow [23, 24]. Hence the result. Theorem 6.2: Let u∗1, u ∗ 2, u ∗ 3 be optimal control functions and S∗h, I ∗ h, R ∗ h,S∗j , I ∗ j ,I∗p are the corresponding state variable of the optimal control problem (2) - (3). Then there exists adjoint variables λ = (λi) T ∈ R6, i = 1, 2, 3, 4, 5, 6, which satisfies the following equations : dλ1 dt = µhλ1 + (1− u1)β1( Ij 1 + kIh )(λ1 − λ2)− (1− u1)β2 SjIh N2 h (λ4 − λ5) dλ2 dt = −X1 + (1− u1) β1ShIjk (1 + kIh)2 (λ2 − λ3) + ηu3γh (1 + αI2h) (λ2 − λ3) + (µh + µ1)λ2 dλ3 dt = δh(λ3 − λ1) + µhλ3 dλ4 dt = −X2 + µjλ4 + (1− u1)β2 Ih Nh (λ4 − λ5) + (1− u2)β3Ip(λ4 − λ5) dλ5 dt = −X2 + (1− u1)β1Sh( 1 1 + kIh )(λ1 − λ2) + µjλ5 dλ6 dt = −X3 + (µp + r − 2rIp k )λ6 (6.7) with transversality conditions λ1(T ) = λ2(T ) = λ3(T ) = λ4(T ) = λ5(T ) = λ6(T ) = 0 (6.8) Proof Let u∗1, u ∗ 2, u ∗ 3 be the optimal control functions and S∗h, I ∗ h, R ∗ h,S∗j , I ∗ j , I ∗ p are the corresponding state variables. Then, Pontryagin’s Maximum Principle ensures the existence of the following adjoint variable λi(i = 1, 2, 3, 4, 5, 6) ∈ R6, which satisfies the following canonical equations: dλ1 dt = − ∂H ∂Sh , dλ2 dt = −∂H ∂Ih , dλ3 dt = − ∂H ∂Rh , dλ4 dt = −∂H ∂Sj , dλ5 dt = −∂H ∂Ij , dλ6 dt = −∂H ∂Ip SENSITIVITY AND OPTIMAL CONTROL ANALYSIS OF JAPANESE 161 with transversality conditions (8), where H is the Hamiltonian defined as above. Thus, the adjoint system (7) can be obtained, In the following result, we shall state the analytical forms of the optimal controls. Theorem 6.3: The optimal controls (u1 ∗, u2 ∗, u3 ∗) which minimizes J over the region U given by u1 ∗ = min{1,max(0, ũ1)} u2 ∗ = min{1,max(0, ũ2)} u3 ∗ = min{1,max(0, ũ3)} where, ũ1 = β1 ( Sh + Ij 1 + kIh ) (λ2 − λ1) + β1ShIjk (1 + kIh)2 (λ3 − λ2) + β2 Ih Nh (λ5 − λ4) Y1 ũ2 = β3Ip(λ5 − λ4) Y2 (6.9) ũ3 = ηγh (1 + αI2h) (λ3 − λ2) Y3 Proof Using the optimal condition (6), we get (9), as follows ∂H ∂u1 = Y1u1 + β1 ( Sh + Ij 1 + kIh ) (λ1 − λ2) + β1ShIjk (1 + kIh)2 (λ2 − λ3) + β2 Ih Nh (λ4 − λ5) = 0 This implies, u1 = β1 ( Sh + Ij 1 + kIh ) (λ2 − λ1) + β1ShIjk (1 + kIh)2 (λ3 − λ2) + β2 Ih Nh (λ5 − λ4) Y1 := ũ1, ∂H ∂u2 = Y2u2 + β3Ip(λ4 − λ5) = 0 This implies, u2 = β3Ip(λ5 − λ4) Y2 := ũ2, And, ∂H ∂u3 = Y3u3 + ηγh (1 + αI2h) (λ2 − λ3) = 0 This implies, u3 = ηγh (1 + αI2h) (λ3 − λ2) Y3 := ũ3. Moreover, lower and upper bounds of these control are 0 and 1 respectively. Thus, if ũ1 > 1, ũ2 > 1, ũ3 > 1, then u1 = u2 = u3 = 1. 162 N.K. GOSWAMI Also, ũ1 < 0, ũ2 < 0 , ũ3 < 0, then u1 = u2 = u3 = 0. Otherwise, we have u1 = ũ1, u2 = ũ2, and u3 = ũ3 Hence, for these controls u∗1, u ∗ 2, u ∗ 3 we get optimum value of the function J . 7. SIMULATION OF OPTIMAL CONTROL In this section, by using MATLAB software, the optimal control model is simulated. The parameter values are remain same as esteemed in the Table 1. The weight constants for the optimal control problem are taken as X1 = 1, X2 = 1, X3 = 1, Y1 = 45, Y4 = 65, Y5 = 75. We solve the optimality system (2) by iterative method with the help of forward and backward difference approximations [23]. We consider the time interval as [0, 200]. Profiles for optimal control u1,u2 and u3 are shown, respectively, in Figures 6(a), 6(b), and 6(c). Finally, to see the effects of optimal controls, the infected human and infected vectors are plotted against time with and without optimal control in Figures 7(b). It is easy to notice that optimal control is more effective in reducing the number of infectives. We consider different types of strategies to see the impact of optimal control in the total number of human infectives. For this purpose, we shall define the following control strategies: • Strategy A: Employing control intervention (u1), only. • Strategy B: Employing control intervention (u2) only. • Strategy C: Employing control intervention (u2) only. • Strategy D: Employing control intervention (u1, u2), (u1, u3), (u2, u3) • Strategy E: Employing control intervention (u1, u2, u3) 7.1. Strategy A: Employing electronic devices, insecticide-treated bed nets, and mosquito repulsive lotions (u1), only. Here, only control measure u1(t) is used to optimize the objective function J , while control intervention u2(t) = u3(t) = 0, were not employed. The influence of u1(t) is demonstrated in Figure 6(a), to minimize the objective function, the optimal control u1(t) is maintained at the maximum level. A single preventive measure can influence the spread of Japanese Encephalitis in the population. 7.2. Strategy B: Increase pig vaccination to reduce the mosquito infection (u2) only. Here, only control measure u2(t) is used to optimize the objective function J , while control intervention u1(t) = u3(t) = 0, were not employed. In Figure 6(b), we present the plots of population and the due to effects of pig vaccination are demonstrated to minimizing the cost and reducing the number of Japanese Encephalitis infections in the population. 7.3. Strategy C: Use of effective medicines as treatment measure of infectious humans (u3) only. Here, only control measure u3(t) is used to optimize the objective function J , while control intervention u2(t) = u3(t) = 0, were not employed. The influence of u3(t) is demonstrated in Figure 6(c), to minimize the disease effect in human population with minimum cost. 7.4. Strategy D: Effective use of control intervention (u1, u2), (u1, u3), (u2, u3). Here, double control interventions (u1(t), u2(t)) when u3(t) = 0, (u1(t), u3(t)) when u2(t) = 0, and (u2(t), u3(t)) when u1(t) = 0, are used to optimize the objective function J. SENSITIVITY AND OPTIMAL CONTROL ANALYSIS OF JAPANESE 163 The control intervention are demonstrated in Figure 6(d), 6(e) and 6(f). It is observed that combination of two control intervention is most effective then using signal control intervention. Significant number of Japanese Encephalitis infective cases reduce. 7.5. Strategy E: Employing all three control interventions (u1, u2, u3). Here all three control interventions (u1(t), u2(t), u3(t)) are used to optimize the objective function J . From Figure 7(a), it is easy to say that by combining all three optimal controls u1(t), u2(t), and u3(t), the total number of infectious individuals decreases significantly. The three optimal control application is the best control strategy to minimize the number of infectives, and will definitely reduce the spread of JEV. 8. CONCLUSION In this study, a non-linear mathematical model is formulated and analyzed for Japanese encephalitis disease with saturated incident ( β1ShIj 1 + kIh ) and saturated treatment function( γhIh 1 + αIh ) . The basic reproduction number is computed by next generation matrix method. Here, it is clearly identify that the parameter β1 is the most sensitive parameter which transmitted Japanese encephalitis disease in human. The sensitivity analysis of the basic reproduction and Japanese encephalitis infected cases are executed in order to determine the relative importance of the model parameters to the disease prevalence. The results of the numerical simulation is executed to support our mathematical findings and demonstrated graphically. Also, we perform here trend of Japanese encephalitis cases of Assam and Uttar Pradesh through line-graph. The infection rate is high in Assam and Uttar Pradesh comparatively other part of the country(India). Three time-department optimal control is incorporated in model system to eliminate the virus from the tropical region using Pontryagin’s maximum principal. Adapting all three optimal control parameters are the best control strategy to minimize the number of infectives, which will reduce the rate transmission of the disease. It is noticed easily that use optimal control in the mode is more effective than without optimal control in reducing the number of infectives in the considered period. REFERENCES 1. World Health Organization, Japanese encephalitis fact sheets. https://www.who.int/news- room/factsheets/ details/Japanese-encephalitis. Accessed 9 May 2019 2. World Health Organization, Vaccine-Preventable Disease. https://www.who.int/immunization/ monitoring-surveillance/burden/vpd/WHO- SurveillanceVaccinePreventable-10-JE-R2. Accessed 9 April 2015 3. Tiwari, S., Singh, R. K., Tiwari, R., and Dhole,T.N.(2012). Japanese encephalitis: a review of the India perspective, The Brazilian Journal of Infectious Diseases.16(6):564- 573. 4. De, A., Maity, K., Jana, S., and Maiti, M.,(2016). Application of various control strategies to Japanese encephalitis: a mathematical study with human, pig and mosquito, Math. Biosci. 282, 46–60. 5. Ndairou, F., Area, I., Torres, D. F. M., (2020). Mathematical Modeling of Japanese encephalitis under Aquatic Environmental Effects. Mathematics. 8 (11), 1880 . 6. Baniya, V., and Keval, R.,(2020). The influence of vaccination on of JE with standard incidence rate of mosquitoes, pigs and humans.Journal of Applied Mathematics and Computing. 64, 519-550 . 164 N.K. GOSWAMI 0 20 40 60 80 100 120 140 160 180 200 Time(in Days) 0 0.2 0.4 0.6 0.8 1 C o n t r o l P r o f il e ( u 1 ) u 1 , when u 2 =u 3 =0 (a) 0 20 40 60 80 100 120 140 160 180 200 Time(in Days) 0 0.1 0.2 0.3 0.4 0.5 0.6 C o n t r o l P r o f il e ( u 2 ) u 2 , when u 1 =u 3 =0 (b) 0 20 40 60 80 100 120 140 160 180 200 Time(in Days) 0 0.2 0.4 0.6 0.8 1 C o n t r o l P r o f il e ( u 2 ) u 3 , when u 1 =u 2 =0 (c) 0 20 40 60 80 100 120 140 160 180 200 Time(in Days) 0 0.2 0.4 0.6 0.8 1 C o n t r o l P r o f il e ( u 1 ,u 2 ) u 1 , when u 3 =0 u 2 , when u 3 =0 (d) 0 20 40 60 80 100 120 140 160 180 200 Time(in Days) 0 0.2 0.4 0.6 0.8 1 C o n t r o l P r o f il e ( u 1 ,u 3 ) u 1 , when u 2 =0 u 3 , when u 2 =0 (e) 0 20 40 60 80 100 120 140 160 180 200 Time(in Days) 0 0.2 0.4 0.6 0.8 1 C o n t r o l P r o f il e ( u 2 ,u 3 ) u 2 , when u 1 =0 u 3 , when u 1 =0 (f) Fig. 7.6. (a) Control profile of u1, when u2 = u3 = 0, (b) Control profile of u2, when u1 = u3 = 0, (c) Control profile of u3, when u2 = u3 = 0, (d) Control profile of u1, u2, when u3 = 0,(e) Control profile of u1, u3, when u2 = 0, (f) Control profile of u2, u3, when u1 = 0 SENSITIVITY AND OPTIMAL CONTROL ANALYSIS OF JAPANESE 165 0 20 40 60 80 100 120 140 160 180 200 Time(in Days) 0 0.2 0.4 0.6 0.8 1 C o n t r o l P r o f il e ( u 1 ,u 2 ,u 3 ) u 1 u 2 u 3 (a) 0 20 40 60 80 100 120 140 160 180 200 Time(in Day) 0 50 100 150 In fe c te d P o p u la ti o n s Without optimal control With optimal control (b) Fig. 7.7. (a) Control profile of u1, u2, and u3, (b) Variation of infected human populations against time with and without optimal control. 7. Panja, P., Mondal, S.K., and Chattopadhyay, J.,(2016).Stability and bifurcation analysis of Japanese encephalitis model with/without effect of some control parameters. Comput. Appl. Math 32(2), 155-175 . 8. Naresh, R. and Panday, S.,(2009). Modelling and analysis of the spread of Japanese encephalitis with environmental effects, Appl. Math. 4(1), 155–175 . 9. Mukhopadhyay, B.B., and Tapaswe, P.K.,(1994). An SIRS epidemic model of Japanese encephalitis. Math. Sci. 17(2), 347–356. 10. Srivastav, A. K., and Ghosh, M.,(2021). Modeling the transmission dynamics of malaria with saturated treatment : a case study of India. Journal of Applied Mathematics and Computing.https://doi.org/10.1007/s12190-020-01469-7. 11. Goswami, N.K.,(2021).Modeling analysis of Zika virus with saturated incidence using optimal control theory. International journal of dynamical systems and differential equations, 11(3-4), 287-301. 12. Gumel, A.B., (2012), Causes of backward bifurcations in some epidemiological models, J.Math. Anal. Appl. 395, 355-365. 13. Chitnis, N., Hyman, J. M., and Cushing, J.M. . Determining important parameters in the spread of malaria through the sensitivity analysis of a mathematical model, Bull.Math.Biol. 70(5), 12721296. 14. Zhang, J., Jai, J., and Song, X., (2O14). Analysis of an SEIR Epidemic Model with Saturated Incidence and Saturated Treatment Function. The Scientific World Journal. vol.2014, ID 910421, 11 Pages. https: //doi.org/10.1155/2014/910421. 15. Goswami, N.K., and Shanmukha, B.,(2020) Stability and optimal control analysis of zika virus with saturated incidence rate. Malaya Journal of Mathematik. 8(2).331-342. 16. Dubey, B., Patra, A., Srivastava, P.K. and Dubey, U. S.,(2013) Modeling and analysis os an SEIR model with different types of nonlinear treatment rates, J. Biol. Syst. 21(3), 1350023. https://doi.org/10.1142/S021833901350023. 17. Zhang, Z., and Suo, Y., (2010). Qualitative analysis of a SIR epidemic model with saturated treatment rate. J. Appl. Math. Comput. 34,177–194 . 18. Li, X. Z., Li, W. S., Ghosh, M., (2009). Stability and bifurcation of an SIS epidemic model with treatment. Chaos Solit. Fract. 42 (5), 2822–2832. 19. Castillo-Chavez, C., Feng, Z.,and Huang, W., (2002). On the computation of R0 and its role on global stability, Mathematical Approaches for For Emerging and Reemerging Infectious Diseases.Springer-Verlag. 9780387953540, 229-250 . 166 N.K. GOSWAMI 20. Driessche, P.V., and Watmough, J., (2002). Reproduction numbers and sub-threshold endemic equilibria for compartmental models of disease transmission. Mathematical Biosciences. 180, 29-48. 21. Goswami, N. K., and Shanmukha, B., (2020).Modeling and Analysis of Symptomatic and Asymptomatic Infections of Zika Virus Disease with Non-Monotonic Incidence Rate. Appl. Math. Inf. Sci. 14 (4), 655-671. 22. Muniaraj, M. and Rajamannar, V. (2019). Impact of SA 14-14-2 vaccination on the occurrence of Japanese encephalitis in India. Hum vaccin Immunother. 15(4),834-840. doi:10.1080/21645515.2018.1564435. 23. Lenhart, S., and Workman, J.T., (2007)Optimal control applied to biological models, CRC Press 24. Pontryagin, L. S., Boltyanskii, V.G., Gamkrelidze, R.V. and Mishchenko, E.F. (1962) The mathematical theory of optimal processes, Inter science Publishers, Geneva. Introduction The Model Formulation Positivity and Boundedness Basic Reproduction Number Data scenario Sensitivity Analysis of R0 Optimal Control Problem Cost Construction and Characterization of Optima Controls Total Cost Determination Existence and Characterization of Optima Controls Simulation of Optimal Control Strategy A: Employing electronic devices, insecticide-treated bed nets, and mosquito repulsive lotions (u1), only. Strategy B: Increase pig vaccination to reduce the mosquito infection (u2) only. Strategy C: Use of effective medicines as treatment measure of infectious humans (u3) only. Strategy D: Effective use of control intervention (u1, u2),(u1, u3),(u2, u3). Strategy E: Employing all three control interventions (u1,u2 ,u3). Conclusion