INTERNATIONAL JOURNAL OF COMPUTING SCIENCE AND APPLIED MATHEMATICS.
VOL.
NO.
AUGUST 2024
Design of Monkeypox Virus Spread Control in Humans Using Pontryagin Minimum Principle Lukman HanafiO .
Mardlijah.
Daryono Budi Utomo.
Suhud Wahyudi and Alya Nur Sha-brina AbstractAiMonkeypox is a contagious disease caused by a In Africa, monkeypox results in death in 1 out of 10 infected individuals.
The Food and Drug Administration in the United States recommends vaccination as a preventive measure against monkeypox virus.
If infected, the World Health Organization (WHO) advises quarantine to prevent further transmission to others.
This research develops a mathematical model known as SIR (Susceptible-Infected-Recovere.
for the spread of monkeypox virus, incorporating vaccination and quarantine as control The SIR model utilized is based on an existing model and follows the conditions of monkeypox spread in Nigeria, represented as a system of nonlinear differential equations.
Optimal control is determined using the Pontryagin Minimum Principle and simulated using the fourth-order forward-backward sweep Runge-Kutta method to assess the level of monkeypox infection before and after implementing control measures.
Based on the simulation results, it is concluded that the application of control measures can reduce the population of infected monkeys by 70% and infected humans by 59%.
INTRODUCTION
INCE January 2022, 3413 laboratories from 50 countries have reported the emergence and fatalities associated with monkeypox, making it a matter of special concern .
Monkeypox is a contagious disease caused by a virus.
The virus is typically transmitted to humans through contact with infected pets or primates, consumption of infected animal meat, direct contact, or animal scratches or bites.
Human-tohuman transmission mainly occurs through direct contact.
efforts to prevent the spread of monkeypox virus, the U.
Food and Drug Administration recommends the use of the JYNNEOSTM vaccine, which has an effectiveness rate of For individuals already infected, quarantine is essential to break the chain of transmission.
Additionally, maintaining personal hygiene, wearing masks, washing utensils with hot water, and disinfecting contaminated surfaces are advised .
This research optimizes the impact of vaccination and quarantine on the spread of monkeypox virus using the Pontryagin Maximum/Minimum Principle.
The Pontryagin Maximum/Minimum Principle is a principle used to solve optimal control problems in the SIR model by finding controls that maximize or minimize the objective function.
It has the advantage of stating the necessary conditions to obtain the most optimal control, thereby minimizing the objective function .
, .
The research also utilizes the Runge-Kutta method to solve the problem numerically.
The Runge-Kutta Hanafi.
Mardlijah.
Utomo.
Wahyudi.
Sha-brina are with the Mathematics Department Institut Teknologi Sepuluh Nopember Surabaya.
Indonesia e-mail: lukman@matematika.
Manuscript received July 18, 2023.
accepted September 14, 2023.
method is an alternative method to Taylor series that does not require derivative calculations .
, .
Its advantages include higher accuracy compared to EulerAos method.
HeunAos method, and Taylor series .
Theoretical studies on monkeypox are conducted by modeling them as systems of differential equations.
Bhunu et al.
conducted research by modeling the spread of diseases like The study concluded that the high transmission rate of monkeypox virus in Central and West Africa is attributed to poor nutrition and poverty, which forces people to hunt monkeys and wild animals that are infected with The study also suggests that further research should estimate the impact of vaccination in reducing monkeypox transmission .
Another study titled AuThe Transmission Potential of Monkeypox Virus in Human PopulationsAy by Fine et al.
concluded that without appropriate interventions, monkeypox has the potential to become a global health threat .
Considering the numerous reports from laboratories regarding monkeypox cases, it is crucial to control the spread of this disease.
Based on the recommendations from BhunuAos research and the preventive measures suggested by the WHO, this study will apply the Pontryagin Minimum Principle to suppress the spread of monkeypox virus.
It is hoped that by implementing vaccination and quarantine controls, the spread of monkeypox virus can be reduced, particularly among susceptible populations.
II.
M ODEL F ORMULATION
Description of the Model This research develops a model for the spread of monkeypox virus, which has been previously studied by Bhunu et al.
, as shown in the following equations:
= un Oe (AAn n )Sn = n Sn Oe (AAn dn An )In = An In Oe AAn Rn = uh Oe (AAh h )Sh = h Sh Oe (AAh dh Ah )Ih = Ah Ih Oe AAh Rh Furthermore, in this study, additional control variables are introduced, including the vaccination rate among susceptible INTERNATIONAL JOURNAL OF COMPUTING SCIENCE AND APPLIED MATHEMATICS.
VOL.
NO.
AUGUST 2024
1 ), the proportion of infected monkeys placed under quarantine .
2 ), and the vaccination rate among susceptible monkeys .
3 ).
As a result, the model is formulated as follows:
= uh Oe .
Ih .
Oe u2 )n2 In u1 AAh ] Sh = .
Ih .
Oe u2 )n2 In ] Sh Oe (AAh dh Ah )Ih = Ah Ih u1 Sh Oe AAh Rh = un Oe [.
Oe u2 )n1 In u3 AAn ] Sn = [.
Oe u2 )n1 In ] Sn Oe (AAn dn An )In = u3 Sn An In Oe AAn Rn Where:
A uh : Birth rate of humans A un : Birth rate of monkeys A AAh : Death rate of humans A AAn : Death rate of monkeys A dh : Death rate due to monkeypox in humans A dn : Death rate due to monkeypox in monkeys A Ah : Natural recovery rate of humans A An : Natural recovery rate of monkeys A n1 : Transmission rate of monkeypox in monkeys from contact with monkeys A n2 : Transmission rate of monkeypox in humans from contact with monkeys A h : Transmission rate of monkeypox in humans from contact with humans A Sh : Population of susceptible humans A Ih : Population of infected humans A Rh : Population of recovered humans A Sn : Population of susceptible monkeys A In : Population of infected monkeys A Rn : Population of recovered monkeys The Equilibrium Points The equilibrium points of the system of equations are obtained when:
The disease-free equilibrium point is a condition where there is no spread of monkeypox virus in a population, resulting in no infected population (Ih = 0.
In = .
Thus, the disease-free equilibrium point can be obtained as follows:
E0 =
u1 Ov uh
u3 Ok n un
, 0,
, 0,
u1 Ov AAh AAh .
1 Ov AAh ) u3 Ok AAn .
3 Ok n AAn )AAn The endemic equilibrium points are used to indicate the potential occurrence of disease transmission.
In essence, there are three possible endemic equilibrium states mathematically:
the specific monkey endemic equilibrium, the specific human endemic equilibrium, and the equilibrium state where the disease coexists between humans and animals.
Based on the fact that monkeypox infection is primarily transmitted from animals to humans, analyzing the endemic equilibrium solely in humans is not necessary since human-to-human transmission of monkeypox rarely causes outbreaks.
The specific animal endemic equilibrium point occurs when there is only infection from animal to animal, no infection from human to human, and no infection from animal to human .
2 = h = .
Therefore, the endemic equilibrium point specific to animal disease can be obtained as follows:
E1O = (ShO , 0.
ROh .
SnO .
InO .
ROn ) uh ShO = u1 AAh u1 ShO ROh = AAh un Oe (An AAn dn )InO SnO = u3 AAn InO = O 3 n An In ROn = AAn The endemic equilibrium points in humans and monkeys occur when there is transmission from animal to animal, animal to human, and human to human.
The endemic equilibrium points in this model can be obtained as follows:
E2O = (ShO .
IhO .
ROh .
SnO .
InO .
ROn ) uh Oe (AAh dh Ah )IhO ShO = u1 AAh uh Oe AAh Ih = O h u1 Sh ROh = AAh un Oe (An AAn dn )InO O Sn = u3 AAn InO = O 3 n An In ROn = AAn The Basic Reproduction Number The basic reproduction number, denoted as R0 , is the expected number of infections generated by a single infected individual in a susceptible population within a unit of time.
In this research, the analysis of the reproduction number is conducted without implementing any control measures.
The determination of the basic reproduction number is performed INTERNATIONAL JOURNAL OF COMPUTING SCIENCE AND APPLIED MATHEMATICS.
VOL.
NO.
AUGUST 2024
using the Next Generation Matrix (NGM) method.
The Next Generation Matrix is defined as follows:
u u K = FV Oe1 = h n1 AAh .
h AAh Ah ) AAh .
h AAh Ah ) un n1 AAn .
n AAn An ) Therefore, the basic reproduction number of the monkeypox disease model is given by:
R0 = {R0,n .
R0,h } where R0,n and R0,h are the basic reproduction numbers for monkeys and humans, respectively, with the following values:
un n1 R0,n = AAn .
n AAn An ) uh h R0,h = AAh .
h AAh Ah ) The disease-free equilibrium point will be asymptotically stable if R0 < 1 and unstable if R0 > 1, where R0 = max{R0,n .
R0,h }.
Before conducting the control analysis, the model is constructed by adding the desired controls.
In this study, controls u1 and u2 are given, which represent vaccination and quarantine controls, respectively.
The controls are restricted to 0 O u O 1.
This leads to the state-space representation, and the matrices A and B are obtained as follows.
Oea Oea Oea OeAA Oea a1 a2 a4 OeAAh OeAh Oedh 0 Ah OeAAh 0 Oea5 OeAAn Oea8 0 Oea6 a6 OeAAn OeAn Oedn Where:
a1 = h Ih a2 = .
Oe u2 )n2 In a3 = u1 a4 = h Sh a5 = .
Oe u2 )n1 In a6 = .
Oe u2 )n2 Sh a7 = .
Oe u2 )n1 Sn Local Stability a8 = u3 Local stability refers to the stability of a linear system or the stability of the linearization of a nonlinear system.
The local stability at an equilibrium point is determined by the signs of the real parts of the characteristic roots of the systemAos Jacobian matrix computed around the equilibrium point.
the case of a nonlinear system, it needs to be linearized to obtain a linear system form.
The following are the properties of local stability around an equilibrium point.
An equilibrium point is said to be asymptotically stable if and only if Re.
) < 0 for every i = 1, 2, .
, n.
An equilibrium point is said to be stable if and only if Re.
) O 0 for every i = 1, 2, .
, n.
An equilibrium point is said to be unstable if and only if Re.
) > 0 for every i = 1, 2, .
, n.
To analyze the stability, we will determine the equilibrium points of the dengue fever transmission model.
The equilibrium points obtained are the disease-free equilibrium E0O = (ShO , 0.
ROh .
SnO , 0.
ROn ), the endemic equilibrium specific to animals E1O = (ShO , 0, 0.
SnO .
InO .
ROn ), and the endemic equilibrium for both animals and humans E2O = (ShO .
IhO .
ROh .
SnO .
InO .
ROn ).
Next, stability analysis is performed by finding the eigenvalues around the equilibrium points.
It is found that the system is asymptotically stable towards the disease-free equilibrium and the endemic equilibrium if all eigenvalues have negative real parts ( < .
OeSh In Sh n2 OeIn Sh n2 In Sn n1 OeIn Sn n1 Next, the rank of Mc = [B | AB | A2 B | A3 B | A4 B | A5 B] is calculated to analyze controllability.
If rank(Mc ) = 6, then the monkeypox spread model is a controllable system.
R ESEARCH M ETHOD
The research was conducted using the following steps:
Literature review on mathematical modeling of monkeypox virus spread.
Modification of the mathematical model of monkeypox virus with vaccination and quarantine.
Determination of equilibrium points and analysis of their stability in the modified mathematical model of monkeypox virus spread.
Control solution using the PontryaginAos Minimum Principle (PMP) method.
Numerical simulation of the modified mathematical model of monkeypox virus spread.
Drawing conclusions and preparing the final report.
Controllability Analysis A system can be controlled if, based on control analysis, it is deemed controllable.
The necessary and sufficient condition for a controllable system is as follows:
,t1 ) = Z t1 eOeAT BBT eOeA T dT is non-singular The matrix Mc = (B | AB | A2 B | A A A | AnOe1 B) has the same rank as n.
IV.
R ESULTS AND D ISCUSSIONS
Formulation and Solution of the Optimal Control Problem The formulation of an optimal control problem consists of mathematically describing a system or model, determining an objective function, and specifying constraints or boundary conditions, with the aim of finding the value of u.
that can optimize the objective function.
INTERNATIONAL JOURNAL OF COMPUTING SCIENCE AND APPLIED MATHEMATICS.
VOL.
NO.
AUGUST 2024
Sn .
n Oe (.
Oe uO2 )n1 In uO3 AAn ) Sn ] In this study, the objective function to be minimized is:
1 , u2 ) = Z tf In [.
Oe uO2 )n1 In Sn Oe (AAn dn An ) In ] Rn .
O3 Sn An In Oe AAn Rn ] Ih .
In .
C1 u21 .
C2 u22 .
C3 u23 .
where Ih represents the population of infected humans and In represents the population of infected monkeys.
C1 is the cost weight coefficient for human vaccination control.
C2 is the cost weight coefficient for infected monkey quarantine control, and C3 is the cost weight coefficient for vulnerable monkey quarantine control.
In other words, the infected population will be minimized by implementing vaccination and quarantine measures with minimum cost.
The Minimum PontryaginAos Principle is used to obtain the optimal control in a dynamic system from the initial state to the final state by minimizing the objective function with control u.
restricted to u.
OO U.
The steps to solve the optimal control problem using the PontryaginAos Minimum Principle are as follows.
Step 1: Formulate the Hamiltonian function H = Ih In C1 u21 C2 u22 C3 u23 Sh .
h Oe .
Ih .
Oe u2 )n2 In u1 AAh ) Sh ] Ih [.
Ih .
Oe u2 )n2 In ) Sh Oe (AAh dh Ah ) Ih ] Rh [Ah Ih u1 Sh Oe AAh Rh ] Sn .
n Oe (.
Oe u2 )n1 In u3 AAn ) Sn ] In [.
Oe u2 )n1 In Sn Oe (AAn dn An ) In ] Rn .
3 Sn An In Oe AAn Rn ] Step 2: Minimize H to all control vectors u.
to determine the stationary conditions.
OCH
OC u1 C1 u1 (Rh Oe Sh )Ov Sh = 0 (Sh Oe Rh )Sh OCH OC u2 C2 u2 (Sh Oe Ih ).
2 In Sh ) (Sn Oe In ).
1 In Sn ) = 0
(In Oe Sn )n1 In Sn (Sh Oe Ih )n2 In Sh
OCH
OC u3 C3 u3 (Rn Oe Sn )Sn = 0 (Sn Oe Rn )Sn Step 3: Use the result from Step 2 by substituting it into Step 1 and determine the optimal H .
H O = Ih In C1 uO2 1 C2 u2 C3 u3 Sh .
h Oe .
Ih .
Oe uO2 )n2 In uO1 AAh ) Sh ] Ih [.
Ih .
Oe uO2 )n2 In ) Sh Oe (AAh dh Ah ) Ih ] Rh [Ah Ih uO1 Sh Oe AAh Rh ] Step 4: Solve the state equations OCH Soh = = uh Oe .
Ih .
Oe u2 )n2 In u1 AAh ) Sh .
OC Sh
OCH
= .
Ih .
Oe u2 )n2 In ) Sh Oe (AAh dh Ah ) Ih .
Ioh = OC Ih
OCH
Roh = = Ah Ih u1 Ov Sh Oe AAh Rh .
OC Rh
OCH
Son = = un Oe (.
Oe u2 )n1 In u3 AAn ) Sn .
OC Sn
OCH
Ion = = .
Oe u2 )n1 In Sn Oe (AAn dn An ) In .
OC In OCH = An In u3 Sn Oe AAn Rn .
Ron = OC Rn And the costate OCH = (Sh Oe Ih ) .
Ih .
Oe u2 )n2 In ] (Sh Oe Rh )u1 AAh Sh .
OC Sh
OCH
IAh = Oe = Oe1 (Sh Oe Ih ) .
Sh ] (AAh dh Ah )Ih Oe Ah Rh .
OC Ih
OCH
RA h = Oe = AAh Rh .
OC Rh
OCH
SA n = Oe = (Sn Oe In ) [.
Oe u2 )n1 In ] (Sn Oe Rn )u3 Sn AAn .
OC Sn
OCH
IAn = Oe = Oe1 (Sh Oe Ih ) [.
Oe u2 )n2 Sh ] (Sn Oe In ) [.
Oe u2 )n1 Sn ] OC In (AAn dn An )In Oe An Rn .
OCH
Rn = Oe = AAn Rn .
OC Rn SA h = Oe Numerical Solution The Runge-Kutta method is a numerical method used to solve initial value problems in differential equations.
The Runge-Kutta method provides smaller errors compared to other numerical methods such as the Euler method and the Heun method.
The fourth-order Runge-Kutta method is widely used because it offers higher accuracy.
LetAos consider the following differential equation as an example:
= f .
, .
In the fourth-order Runge-Kutta method, it is formulated as yn 1 = yn .
1 2k2 2k3 k4 ) k1 = f .
n , yn ), k2 = f xn , yn k3 = f xn , yn INTERNATIONAL JOURNAL OF COMPUTING SCIENCE AND APPLIED MATHEMATICS.
VOL.
NO.
AUGUST 2024
k4 = f .
n h, yn k3 ).
To obtain the optimal control, a numerical solution using the fourth-order forward-backward sweep Runge-Kutta method is This is done because the state equation is given by Sh .
= Sh0 .
Sn .
= Sn0 .
Ih .
= Ih0 .
In .
= In0 .
Rh .
= Rh0 .
Rn .
= Rn0 , while the costate equation is given by the final values 1 (T ) = 2 (T ) = 3 (T ) = 4 (T ) = 5 (T ) = 6 (T ) = 0.
Analysis of Simulation Results In this discussion, the initial conditions for each population are given as follows:
Sh .
= 0.
Ih .
= 0.
Rh .
= 0.
Sn .
= 0.
In .
= 0.
Rn .
= 0.
According to the development of monkeypox cases in Nigeria, the parameters of this study are obtained from the research conducted by Bhunu et al.
with the following Fig.
2: Graph of Infected Human Population that implementing control measures reduces the population by 40%.
This decrease occurs because the control measures inhibit the spread of the virus.
TABLE I: Model Parameter Values Parameter un uh AAn AAh An Ah n1 monkeys from contact with monkeys n2 humans from contact with monkeys h humans from contact with humans Value 2 yrOe1 029 yrOe1 15 yrOe1 02 yrOe1 3 yrOe1 33 yrOe1 2 yrOe1 1 yrOe1 87 yrOe1 Description birth rate of monkeys birth rate of humans death rate of monkeys death rate of humans natural recovery rate of monkeys natural recovery rate of humans death rate due to monkeypox in monkeys death rate due to monkeypox in humans .
transmission rate of monkeypox in 62 yrOe1 .
transmission rate of monkeypox in 73 yrOe1 .
transmission rate of monkeypox in With the simulation results as follows:
Fig.
1: Graph of Infected Monkey Population Figure 1 and 2 show a decrease in the population of infected humans and monkeys.
The simulation results indicate Fig.
3: Graph of Susceptible Monkey Population From Figure 3, the monkey population experiences a decline in the first two years due to the movement of vulnerable monkey population towards the recovered monkey population through quarantine measures.
However, the population starts to increase as the control over quarantine measures decreases.
From Figure 4, there is a decrease in the susceptible human population from the first year to the fourth year.
This decline in population in the simulation is due to the movement of susceptible humans towards the recovered human population through vaccination.
However, after the fourth year, the susceptible population starts to increase as the vaccination control in humans decreases.
From Figure 5, it can be observed that there is an increase in the population of recovered monkeys.
The simulation shows a 20% increase in the population, which is attributed to the influence of the movement from the vulnerable monkey population towards the recovered monkey population through quarantine measures.
INTERNATIONAL JOURNAL OF COMPUTING SCIENCE AND APPLIED MATHEMATICS.
VOL.
NO.
AUGUST 2024
In Figure 6, there is a 5% increase in the population in the This increase is due to the movement from the susceptible human population towards the recovered human population through vaccination.
Fig.
4: Graph of Susceptible Human Population Fig.
7: Graph of Recovered Monkey Population Figure 7 shows that the green line represents the vaccination rate .
1 ), the blue line represents the quarantine control on infected monkeys .
2 ), and the yellow line represents the quarantine control on vulnerable monkeys .
3 ).
In the year when vaccination and quarantine control measures are implemented, the required vaccination rate for the susceptible human population is 1 per year, the proportion of infected monkeys to be quarantined is 0.
1, and the required quarantine rate for the vulnerable monkey population is 0.
66 per year.
C ONCLUSIONS
Fig.
5: Graph of Recovered Monkey Population Based on the previous analysis and discussion, several conclusions can be drawn as follows:
The mathematical model developed for the simulation is as follows:
= uh Oe .
Ih .
Oe u2 )n2 In u1 AAh ] Sh , = .
Ih .
Oe u2 )n2 In ] Sh Oe (AAh dh Ah ) Ih , = Ah Ih u1 Sh Oe AAh Rh , = un Oe [.
Oe u2 )n1 In u3 AAn ] Sn , = [.
Oe u2 )n1 In ] Sn Oe (AAn dn An ) In , = u3 Sn An In Oe AAn Rn .
By using the Maximum PontryaginAos Principle, the optimal control u obtained from the mathematical model of monkeypox spread is as follows:
Fig.
6: Graph of Recovered Human Population (Sh Oe Rh )Sh (In Oe Sn ).
1 In Sn ) (Ih Oe Sh ).
2 In Sh ) INTERNATIONAL JOURNAL OF COMPUTING SCIENCE AND APPLIED MATHEMATICS.
VOL.
NO.
AUGUST 2024
(Sn Oe Rn )Sn .
The simulation results after implementing the control
measures indicate that the population of infected monkeys and humans can be reduced by 40
ACKNOWLEDGMENT