Modelling the dynamics of Pine Wilt Disease with asymptomatic carriers and optimal control

Pine wilt disease is a lethal tree disease caused by nematodes carried by pine sawyer beetles. Once affected, the trees are destroyed within a few months, resulting in significant environmental and economic losses. The role of asymptomatic carrier trees in the disease dynamics remains unclear. We developed a mathematical model to investigate the effect of asymptomatic carriers on the long-term outcome of the disease. We performed a stability and sensitivity analysis to identify key parameters and used optimal control to examine several intervention options. Our model shows that, with the application of suitable controls, the disease can be eliminated in the vector population and all tree populations except for asymptomatic carriers. Of the possible controls (tree injection, elimination of infected trees, insecticide spraying), we determined that elimination of infected trees is crucial. However, if the costs of insecticide spraying increase, it can be supplemented (although not replaced entirely) by tree injection, so long as some spraying is still undertaken.

www.nature.com/scientificreports www.nature.com/scientificreports/ Here, we develop a dynamic model of PWD incorporating an asymptomatic carrier class and examine control policies that minimize implementation costs while protecting forests from the disease. To the best of our knowledge, none of the previous mathematical studies used optimal control to explore the transmission dynamics of the PWD in the presence of the asymptomatic carriers.

Model formulation
The total host (pine wood trees) and vector (beetles) are represented by N H (t) and N V (t), respectively. N H (t) is further classified into four epidemiological classes: susceptible pine trees S H (t), exposed pine trees E H (t), asymptomatic carrier pine trees A H (t) and infected pine trees I H (t). N V (t) is classified into three epidemiological classes: susceptible beetles S V (t), exposed beetles E V (t) and infected beetles I V (t).
The recruitment rates of host trees and beetles are represented, respectively, by Λ H and Λ V , while the natural death rates of host pine trees and vector beetles are denoted by γ 1 and γ 2 , and the disease mortality rate of host pine trees is represented by μ. Here, m and η are the respective rates of progression from the exposed class to the infected class in the host and vector populations. The term β 1 ψS H I V denotes the incidence rate, where β 1 is the rate of transmission and ψ is the average number of daily contacts with vector adult beetles during maturation. β 2 is the rate at which an infected beetle transmits a nematode through oviposition, with the average number of oviposition contacts per day denoted by θ. The termination of oleoresin exudation in susceptible trees without infection of nematode is denoted by α. We thus interpret β 2 θα as the transmission through oviposition, and hence β 2 θαS H I V represents the number of new infections. A fraction ω (0 ≤ ω ≤ 1) of the exposed tree class    www.nature.com/scientificreports www.nature.com/scientificreports/ The total population sizes of host and vector are given by For biological realism, we study the behaviour of the system (1) in the closed set Nonnegative solutions of system (1) can be easily verified for appropriate initial values. The first four equations of (1) imply that  Similarly, adding the last three equations of the system (1), we get Using the comparison theorem again, there exists t 2 > t 1 , such that Hence, the solutions of the system (1) are bounded.
In the Supplementary Material, we determine 0  and prove that the disease-free equilibrium (DFE) is globally asymptotically stable, which also rules out the possibility of a backward bifurcation. We also show that the endemic equilibrium is globally asymptotically stable, under certain conditions.

Sensitivity analysis of threshold dynamic
Due to uncertainties in experimental data, determining accurate outcomes from an epidemiological system is difficult 25 . To compensate for these uncertainties, we use partial rank correlation coefficients (PRCCs) to identify the impact of all parameters on  0 . This technique measures the degree of the relationship between inputs and output of the system. Positive PRCCs indicate parameters that increase  0 when they are increased, while nega- www.nature.com/scientificreports www.nature.com/scientificreports/ tive PRCCs indicate parameters that decrease  0 when they are increased. Parameters with PRCCs values greater than 0.4 in magnitude have a significant effect on the outcome. Figure 2 illustrates the effect of parameter variations on  0 for all fourteen parameters. Clearly, 0  is most sensitive to γ 1 and γ 2 , the natural death rates of pine trees and beetles, respectively; the latter can be controlled using insecticide (u 3 ), while the former can be partially controlled by eliminating infected trees (u 2 ). 0  is also sensitive to the birth rates of pine trees and beetles, the latter of which can be controlled using insecticide (u 3 ). The transmission rate K is also a sensitive parameter, which can be controlled by nematicide-injection and vaccination (u 1 ).

optimal control strategies
In this section, we introduce u 1 , u 2 and u 3 as three control measures that can affect PWD. The force of infection in the pine-tree population is reduced by (1 − u 1 ), where precautionary measures efforts are denoted by u 1 ; for example nematicide injection and vaccination. To keep the host tree population safe and to prevent infection, the nematicide-injection preventative control measure is used. We use the control variable u 2 to describe elimination of infected host trees. Supplementary infections are extremely reduced by demolition and elimination of infected host trees. The removal of these infected trees guarantees that eggs, larvae and pupa that are occupying the host pines are devasted. Our third control variable represents spraying of insecticide and larvacide to kill adult insects and reduce the vector birth rate.
Model (1) is modified for optimal control as follows: ,  www.nature.com/scientificreports www.nature.com/scientificreports/ with nonnegative initial conditions. The control functions u(t) = (u 1 , u 2 , u 3 ) ∈ U associated to the variables S H , E H , The constants b 0 and b 1 are removal-rate constants whose inverses correspond to the average time spent in the relevant compartment. Since it is unlikely that infected trees will be removed within one day of infection, we set b 1 = 1; hence the range 0 ≤ u 2 ≤ 1 corresponds to a removal time between 1 day and infinite time. The objective functional for the optimal-control problem is subject to the control system (2). The constants L 1 , L 2 , L 3 , L 4 , B 1 , B 2 and B 3 are the weight or balancing constants, which measure the relative cost of interventions over the interval [0, T]. We seek optimal controls ⁎ ⁎ ⁎ u u u , , Clearly, the equations in the control system (2) are bounded above, and thus we can apply the results in 26 to model (2). Moreover, the set of control variables and the state variables is nonempty, and the set of control variables denoted by U is closed and convex. In the control problem (2), the right-hand side is continuous and bounded above by state variables and a sum of the bounded control, and can be expressed as a linear function of U having state-and time-dependent coefficients. Hence there exists constants m > 1 and l 1 , l 2 > 0 such that the integrand L(y, u, t) of the objective functional J is convex and satisfies We apply the results presented in 27 to justify the existence of (2) and to obey the above conditions. Clearly, the set of control and state variables are bounded and nonempty. The solutions are bounded and convex. Therefore the system is bilinear in the control variables. We verify the last condition: 3 . In order to get the solution of the control problem, it is necessary to obtain the Lagrangian and the Hamiltonian of (2). The Lagrangian L is expressed as By choosing X = (S H , E H , I H , S H , E H , I H ), U = (u 1 , u 2 , u 3 ) and λ = (λ 1 , λ 2 , λ 3 , λ 4 , λ 5 , λ 6 , λ 7 ), the Hamiltonian can be written   We use Pontryagin's Maximum Principle 28 to obtain the optimal solution of the control system (2). Since   (2), there exist adjoint variables λ i (i = 1, 2, 3, 4, 5, 6, 7) such that www.nature.com/scientificreports www.nature.com/scientificreports/ ( ) ( ) ( ) ( ) ( ) ( ) ( ) Furthermore, the controls ⁎ ⁎ ⁎ u u u , , 1 2 3 are given by  www.nature.com/scientificreports www.nature.com/scientificreports/ Proof. To determine the required adjoint system (8) and the transversality conditions mentioned in (9), we utilize the Hamiltonian in (6). By applying the third condition of (7), we get (8). Applying the second condition of (7), we get (9).

Elimination of infected trees (u 2 ) and spraying of insecticides (u 3 ).
We considered two controls: the elimination of infected trees (u 2 ) and the spraying of insecticides (u 3 ) in the absence of tree injection and vaccination. Figures 3 and 4 show the outcomes in both the absence and presence of control. Figure 3 shows the dynamics of the pine-tree population, while Fig. 4 shows the dynamics of the vector population. With these controls, we see a rapid increase in the population of susceptible trees (Fig. 3(a)) and eventual elimination of exposed and infected trees ( Fig. 3(b,d)), with only the asymptomatic carriers remaining in the infected classes (Fig. 3(c)). The vector population is eventually depleted (Fig. 4(a-c)) in the presence of these two controls. The two control profiles u 2 and u 3 are bounded up to 0.4 and 0.8 (Fig. 4(d)). Biologically, u 2 is the additional elimination rate of only infected trees, while u 3 acts to simultaneously increase the removal rate of all vectors, while also decreasing the birth rate. Since all interventions range between 0 (no control) and 1 (complete control), this suggests that our objective can be achieved with only partial controls. Hence if infected trees are removed 2.5 days or later after infection or if insecticides/larvacides are up to 80% effective, the infection can be controlled.
Tree injection (u 1 ) and spraying of insecticides (u 3 ). We next examine the combination of tree injection (u 1 ) and insecticide spraying (u 3 ). The results are shown in Figs. 5 and 6. With these two controls, there is a significant increase in the population of susceptible and exposed pine trees, while the population of asymptomatic carriers and infected pine trees are reduced but not eliminated (Fig. 5). This suggests that the elimination of infected pine trees has a significant impact on the disease. Note that the vector population is eliminated using these controls (Fig. 6). www.nature.com/scientificreports www.nature.com/scientificreports/ Tree injection (u 1 ) and elimination of infected trees (u 2 ). Considering u 1 and u 2 in combination, Figs. 7 and 8 illustrate that, without insecticide spraying, the control (minimization and/or elimination) of infection in the pine trees is not possible. While the population of susceptible pine trees has a slower decline with these control (Fig. 7(a)), the infection eventually takes over. Likewise, although the susceptible beetle population is recovered using these controls, the infection nevertheless eventually dominates (Fig. 8). It follows that, without insecticide spraying, the control of infection is not possible.
Complete control. We now apply all three controls in order to determine the ideal outcome (Figs. 9 and 10).
Comparing Fig. 9 to Fig. 3, we see that susceptible pine trees recover faster and the disease is eliminated quicker, except for asymptomatic carriers. We thus see that the most effective strategy is to apply all three controls, although similar results can be achieved by applying only two controls: elimination of the infected pine trees (u 2 ) and the spraying of insecticides (u 3 ).
From Fig. 11, we see that, as the cost of u 1 increases, the control profile is dominated by u 3 . That is, if tree injection becomes prohibitively expensive, the procedure can be replaced by increased insecticide spraying. Figure 12 shows little variation in the control profiles as the cost B 2 increases unless the cost is prohibitive. This suggests that the control u 2 is worth implementing, even at high cost. The combination of u 1 and u 3 alone does not eliminate infection, so it follows that elimination of infected trees is essential to disease control. This may hinder disease eradication if the costs of elimination become prohibitively expensive. Figure 13 shows that if the cost of insecticide spraying increases, the control profile is dominated by tree injection. Interestingly, while the combination of u 2 and u 3 produced superior results to the combination of u 1 and u 2 , the latter combination can still produce effective results if supplemented by a small amount of insecticide spraying. www.nature.com/scientificreports www.nature.com/scientificreports/

Discussion
We developed a mathematical model to examine the effect of asymptomatic carriers of Pine Wilt Disease (PWD) on the long-term course of disease. We showed that the disease-free equilibrium was globally asymptotically stable and that the endemic equilibrium was globally asymptotically stable under some conditions. A sensitivity analysis identified key parameters: natural death rates in trees and beetles; birth rates in both trees and beetles; and transmission rates from trees to beetles.
We applied several controls to our system: tree injection, insecticide spraying and elimination of infected trees. These were chosen in conjunction with the most sensitive parameters except for the natural birth and death rates of trees, since our ultimate goal is the preservation of trees. We showed that the disease can be eliminated using suitable controls, except for the asymptomatic carriers. By including this class, our model showed that the disease may remain endemic, requiring permanent control, even in the best-case scenario.
Examining the controls in detail, we found that elimination of infected trees is critical, especially when used in conjunction with insecticide spraying. If the cost of insecticide spraying becomes prohibitive, it can be partially replaced by tree injection. However, if the costs of elimination of infected trees becomes prohibitive, there is no simple replacement, which may result in runaway costs.
It follows that we can control the disease using suitable interventions, but it will not be eliminated due to the presence asymptomatic carriers. The presence of infection in these carriers suggests that infection can restart in nearby healthy trees. It follows that our control measures must be undertaken continually unless such asymptomatic carriers can be identified and removed. This has long-term implications for disease management and economic investment.