OPTIMAL MULTIPLICATIVE CONTROL OF BACTERIAL QUORUM SENSING UNDER EXTERNAL ENZYME IMPACT

. The use of external enzymes provides an alternative way of reducing communication in pathogenic bacteria that may lead to the degradation of their signal and the loss of their pathogeneity. The present study considers an optimal control problem for the semilinear reaction–diﬀusion model of bacterial quorum sensing under the impact of external enzymes. Estimates of the solution of the controlled system are obtained, on the basis of which the solvability of the extremal problem is proved and the necessary optimality conditions of the ﬁrst-order are derived. A numerical algorithm to ﬁnd a solution of the optimal control problem is constructed and implemented. The conducted numerical experiments demonstrate an opportunity to build an eﬀective strategy of the enzymes impact for treatment.

have been proposed in the form of initial-boundary value problems for a system of parabolic partial differential equations (PDEs) [18,19], time-lagging parabolic partial differential equations [20] as well as fractional partial differential equations [21,22]. In particular, the study [19] has reported the numerical simulation results with a focus on AHL and Lactonase concentrations changing under the external addition of substances. In this case, the direct problem has been solved to confirm that the repeated addition of Lactonase enzymes can significantly reduce the concentration of signaling molecules. Since the AHL concentration is associated with the quorum sensing "level", we can make bacteria "fall silent" (i.e. having the quorum sensing system downregulated) during the period of action of the external Lactonase concentration. In this way, the mathematical model of bacterial quorum sensing may be addressed to predicting inter-cell bacterial communication including problems of controlled degrading signaling compounds.
The use of methods of mathematical modeling in bacterial quorum sensing allows us to estimate the required enzyme concentrations and time duration of their exposure. The optimal control approach makes it possible to apply an algorithm of finding the optimal state instead of enumerating various scenarios of enzymes impact. This is especially important, as the time points of exposure may influence the behaviour of the whole system heavily due to the nonlinearity behind. Although optimal control methods for nonlinear reaction-diffusion models are quite developed (see, e.g., [1,5,6,17]), nevertheless, for bacterial quorum sensing model this approach requires a separate study. In spite of wide spectrum of studies devoted to the development of the models of bacterial communication, there is no previous study using optimal control techniques with targets of adjustment of bacterial quorum sensing in the spatially heterogeneous situation.
In the current work, the optimal control problem for the semi-linear reaction-diffusion model of bacterial quorum sensing is examined. The external enzyme concentration is chosen as control. The minimization of the objective functional leads to a reduction of the AHL-concentration under the simultaneous limitation on the total amount of the external enzymes. The solvability of the optimal control problem is proved and the necessary optimality conditions of the first-order are derived. A numerical algorithm to find a solution of the optimal control problem is constructed and implemented. The convergence of the iteration process is analysed. The conducted numerical experiments demonstrate an opportunity to build an effective strategy of the enzymes impact.

Brief biological setups and mathematical formulation of the optimal control problem
As previously introduced, quorum sensing is one of the mechanisms by which the communication of bacteria can be realized. The quorum sensing is the ability of bacteria to communicate and coordinate their behavior due to the secretion of specific molecular signals caused by genes regulation. This phenomenon was discovered by biologists Bonnie Bassler and Everett Peter Greenberg, which elucidated molecular mechanisms underlying quorum sensing and invented new ways for modulating the microbiome for medicine purposes. In order to introduce the main characteristics of the biological system, let us consider a concise description of the bacterial communication process, without pretending to provide a complete presentation of all biological foundations of this multifaceted and complex phenomenon. To be more specific, we will assume that the class of Gram-negative bacteria is subject to consideration, namely P. putida bacterial strain. This bacterial strain is one of the most studied from the point of view of observation of joint dynamic processes as the quorum sensing and quorum quenching due to the production of a special Lactonase enzyme.
The simplified scheme illustrating the quorum-sensing system in P. putida is shown in Figure 1. Bacterial quorum-sensing is realized by means of the generation and distribution of signaling molecules or autoinducers (N-acyl homoserine lactones or AHL in short notation), which allows the bacterial colony to switch the behaviour dependent on the local population density by gene regulation. In addition, the key elements of the gene regulation system are the following: AHL synthases, which are enzymes responsible for the synthesis of autoinducers (proteins of the so-called LuxI family), and transcriptional regulators, which are proteins of the LuxR family and act as receptors for the AHL molecules. Produced signaling molecules at a certain concentration are able to diffuse through the bacterial cell membranes and regulate gene activity, which again results in an increased AHL production due to positive feedback. During the "evolution" of the bacterial population, natural degradation of the signaling substance AHL can occur. To be definite, we assume that the external enzyme concentration is added to the biological system. For instance, these compounds can be associated with Lactonase enzymes which are able to degrade autoinducers of AHL type.
Following the study [21], the considered process is modeled in the domain Ω ⊂ R 3 with the boundary Γ by the following initial-boundary value problem: ∂y/∂t − a∆y + by = f h(y) − uy, x ∈ Ω, t ∈ (0, T ), (2.1) where y is the AHL concentration in mol/l; the unknown function u(t), t ∈ (0, T ) denotes the control function in 1/h, which can be defined as u = γ · Lc, where γ is the degradation rate of AHL by enzyme in l/(mol h) and Lc is the enzyme concentration in mol/l; where x j are fixed points in Ω. The function h(s), s ∈ R is odd and h(s) = α + βs r s 0 + s r , s ≥ 0.
Positive constants a, b, α, β, σ, s 0 , and r > 1 describe the properties of the model and are determined according to [21]. Nonnegative functions g(x), y b (x), x ∈ Γ, and y 0 (x), x ∈ Ω, are specified according to a valid setup of the model due to empirical data. By ∂ n we denote the derivative in the direction of the outward normal n with respect to domain Ω.
The extremum problem is to find the pair { y, u} such that on solutions of initial-boundary value problem (2.1)-(2.3). The function y d and regularization parameter λ > 0 are given. The main goal of the work is to construct an algorithm for solving the posed optimal control problem. A priori estimates for the solution of the controlled system (2.1)-(2.3) are obtained, on the basis of which the solvability of the extremal problem is proved and an optimality system is constructed. The efficiency of the proposed iterative algorithm for finding the optimal control is illustrated by numerical examples.

Formalization of the problem
By L p (0, T ; X) we denote the space of p-integrable on (0, T ) functions assuming values in a Banach space X. Accordingly, by C([0, T ], X) we denote the space of continuous on [0, T ] functions assuming values in a Banach space X. Define the space where v = dv/dt. Notice that W ⊂ C([0, T ]; H) is the continuous embedding, and moreover, the embedding W ⊂ L 2 (0, T ; H) is compact.
Further, we assume that the following conditions hold: Let us define the operator A : V → V and functional g b ∈ V using the following equalities valid for any y, z ∈ V : The bilinear form (Ay, z) determines the inner product in the space V , and the corresponding norm z V = (Az, z) is equivalent to the standard norm of V .
To formulate the optimal control problem, we define the constraint operator F (y, u) :

Solvability of the problem (CP)
Let us obtain preliminary estimates for the solution of the initial-boundary value problem (2.1)-(2.3). Moreover, Here, Proof. Let us obtain a priori estimates for the solution of the problem (3.1). From the equality (y + Ay, y) = (f h(y) − uy + g b , y), taking into account that |h(y)| ≤ α + β, it follows the estimate Therefore, Integrating the last inequality over t and using the Gronwall's inequality to estimate y(t) 2 , we obtain Thus, the nonlocal estimates (4.1) are true. The unique solvability of the parabolic problem (3.1) with Lipschitz nonlinearity h(y) is proved by a standard way. First, we define for small T 1 > 0 the operator S : C([0, is a solution of the following problem on (0, T 1 ): Further, we verify that the operator S is a contraction. Thus, there is a unique fixed point y of the operator S which is a solution of the problem (3.1) on (0, T 1 ). The a priori estimates (4.1) do not depend on length T 1 of the interval of existence of the local solution. Therefore, the local solution can be continued over the entire interval (0, T ). The uniqueness of the solution follows from the Gronwall's inequality. Proof. Let j = inf J on the set u ∈ U , F (y, u) = 0. Then we choose the minimizing sequences u m ∈ U, y m ∈ W , J(y m , u m ) → j, such that From the boundedness of the sequence u m in U , on the base of Lemma 4.1, it follows the estimates Here, by C > 0 we denote the maximal constant bounding the norms and it is independent of m. From these estimates and equalities (4.2), it follows that the sequence {y m } is bounded in L 2 (0, T ; V ). Taking into account the compactness of the embedding W ⊂ L 2 (0, T ; H), and passing, if necessary, to subsequences, we conclude that there is a pair { y, u} ∈ W × U such that u m → u weakly in U, y m → y weakly in L 2 (0, T ; V ), strongly in L 2 (0, T ; H). Note also that ∀z ∈ V The convergence results (4.3)-(4.5) allow us to pass to the limit in (4.2). Therefore, wherein j ≤ J( y, u) ≤ lim J(y m , u m ) = j. Thus, the pair { y, u} is a solution of the problem (CP).

Optimality conditions
To obtain an optimality system, it is sufficient to use the Lagrange principle for smooth-convex extremal problems [11,13]. Let us verify the validity of the key condition that the image of the derivative of the constraint operator F (y, u) with respect to y, where y ∈ W, u ∈ U , coincides with the space L 2 (0, T ; V ) × H. It is this condition that guarantees the nondegeneracy of the optimality conditions. Recall that F (y, u) = {y + Ay − f h(y) + u y − g b , y(0) − y 0 }. Lemma 5.1. Let the conditions (i) hold. Then for any pair {y, u} ∈ W × U the following equality is true: Assuming in (5.2) that z ∈ C ∞ 0 (0, T ; V ), we conclude that, in the sense of distributions on (0, T ), the following equality hold:  wherein λ u = −( y, p).

Algorithm of solving the problem and its implementation
Let us present an algorithm for solving the control problem. Let J(u) = J(y(u), u), where y(u) is, corresponding to the control u ∈ U , a solution of the problem (2.1)-(2.3). According to (5.3), the gradient of the objective functional is determined by J (u) = λu + (y, p).
Here, u ∈ U is the control, p is the adjoint function satisfying (5.5), where y, u are replaced by y, u.
The proposed algorithm for solving the problem (CP) is as follows: Gradient descent algorithm 1: Choose an initial value of the gradient step ε 0 , 2: Choose a number of iterations N , 3: Choose an initial approximation for the control u 0 ∈ U , 4: for k ← 0, 1, 2, . . . , N do:

5:
For a given u k , calculate the state y k , the solution to the problem (2.1)-(2.3).

6:
Calculate the value of the objective functional J(y k , u k ).

7:
Calculate the adjoint state p k from (5.5), where y = y k , u = u k .
The value of the parameter ε k is chosen empirically to provide the convergence of the iterative algorithm. The number of iterations N is chosen sufficient to satisfy the condition J(y k , u k ) − J(y k+1 , u k+1 ) < δ, where δ > 0 determines the accuracy of calculations.
The examples considered below illustrate the convergence of the proposed algorithm for small values, that is important, of the regularization parameter.

Numerical experiments
In order to perform numerical experiments let us consider the bacterial strain Pseudomonas putida IsoF as an object of mathematical modelling. P. putida is a Gram-negative, rod-shaped bacterium, that lives in soils, waters, and plants, particularly, polluted sites [32]. The US National Institutes of Health assigned P. putida as a safety strain. Indeed, certain strains of P. putida are not pathogenic due to lack of genes that determine virulence. P. putida IsoF is often used for experiments and explorations in the context of Quorum sensing due to its versatility and ease of handling.
But some studies have also found out that P. putida can be isolated as opportunistic human pathogens capable of causing nosocomial infections [10]. Acquired diseases are rather rare and are diagnosed in patients with weakened immune systems such as newborns, and cancer patients. Also there is a study saying, that P. putida can provide an exchange platform for resistance genes, and as a consequence, could cause a hospital spread of more virulent microorganisms such as deadly P. aeruginosa [26].
Most infections caused by these bacterial species exhibit resistance to certain antibiotics and its combinations [24]. Recent studies have indicated the quorum sensing phenomenon as well as the production of a special enzyme (Lactonase) that can degrade signaling molecules in P. putida bacteria [9,29]. Hence, P. putida represents a relevant object to perform simulations of optimal control of quorum sensing in bacterial population by the influence of external compounds.
For the numerical experiments we consider the one-dimensional domain Ω = (0, l), where l = 100 µm. Let us also assume that there is one bacterial colony with linear size of 10 µm located at the central position x c = 50 µm of the model domain. The time of observation is 26 hours. We set the model parameters as shown in Table 1 using previously estimated values (see [18]). The values of y d and y 0 are set to zero. Also, an approximation of the zero Dirichlet boundary conditions is used such that g = 10 5 a and y b = 0.
At the first step of the iterative algorithm, the zero approximation of the control is set, that is u 0 = 0. The numerical solution of the problem (2.1)-(2.3) is based on the three-layer implicit finite-difference scheme with the second order temporal and spatial accuracy. Due to the presence of the non-linear reaction terms in (2.1), an iterative procedure is used to obtain the approximation of the solution for each time layer. To solve the problem (5.5) for the adjoin state, we apply the same numerical method, with the exception that the problem is solved for the reverse time direction. Figure 2 shows the time-dependent profiles of the distributions of AHL concentrations calculated under control compared to the case without control u = 0. AHL concentrations are estimated at the central position x c = 50 µm of the model domain. Here, we performed 1000 iterations at different values of regularization parameter λ. The iteration parameter ε k is specified as 0.01/λ. The coordinate distribution of AHL concentrations at the time moment t * = 24 h is given in Figure 3.
From a practical point of view, it makes sense to examine the total effect of external enzyme impact. Figure Fig. 3, λ = 7 × 10 −13 ). In addition, a particular interest is analyzing the L 2 norm of the AHL concentration for different time moments. For example, at the time moment of 24 hours, the AHL concentration decreases by 6.2 times compared with the case of zero control (see Fig. 4, λ = 7 × 10 −13 ).
The time dependencies of the control function u(t) for the different values of the regularization parameter λ are shown in Figure 5. The control function increases rapidly at the beginning of observation period (up to about 1.5 hours), then it drops sharply and passes to constant level. For example, for λ = 7 × 10 −13 , the control at the constant level is about 70 % of the peak value. This means that at the beginning of the therapeutic process it is required a relatively high dose of the enzyme. This observation fits together with the idea that it is more effective, to avoid already the upregulation of the quorum sensing systems which would be expected otherwise after a few hours than to bring it down later from an already high level. Further, it is enough to     provide a certain supporting dose in order to reach an effective decrease in the AHL concentration. In addition to the stabilizing role, the regularization parameter λ is responsible for limiting the total amount of the enzyme.
Note that, since p(x, T ) = 0, the control u(t) is equal to zero at the final time moment. This leads to an increase in the AHL concentration in a certain neighborhood of the final time moment. For this reason, in the discussion of the numerical results presented in Figures 2-4, we considered the behavior of the AHL concentrations for a period of time up to 24 hours, when the amount of introduced enzymes is maintained at a stable level. Figure 6 demonstrates the numerical convergence of the iterative process. The time-dependent profiles of AHL concentrations are calculated at the central position of solution domain x c = 50 µm under variations of the number of iterations: 100, 500, and 1000. The numerical analysis indicates that no further increase in number of iterations is required to provide an acceptable accuracy.

Conclusion
The control techniques for bacterial quorum sensing underlie the approaches to prevent the cell-to-cell communication processes for many bacterial strains. The potential advantage of such approach consists in application of natural enzymes as an alternatives to antibacterial drugs in the case of resistant pathogens. Obviously, this treatment does not directly address the presence of the bacteria, but could control very well their pathogeneity which is often upregulated by a Quorum sensing system. Nevertheless, it can indirectly limit the bacterial growth by avoiding their quorum sensing upregulation and by that a better adaptation to the given environmental conditions.
The present study provides a theoretical background with respect to using external enzyme compounds to reduce the concentration of signaling molecules responsible for the quorum sensing level of bacterial communities and, as a consequence, for their virulence. For that it is not necessary to bring the autoinducer concentration to zero, but just below a certain threshold level, which is much more realistic.
The antivirulence strategy of bacterial quorum quenching by the addition of external enzymes has been formalized as an optimal control problem. The problem is to minimize the objective functional that allows us to reduce the AHL concentration and simultaneously to limit the total amount of used enzymes. The proposed approach leads to an effective scenario of the enzymes impact, which is characterized by the use of an increased amount of enzymes at the beginning of the treatment and their reduced dosage in the subsequent period. Furthermore, the optimised treatment schedule, which can be determined by the optimal control strategy as introduced above, allows to restrict the total amount of enzyme to a chosen level or total amount, instead of using constant high levels thereof which may lead to undesired side effects. One should keep in mind that also beneficial bacteria may use similar kinds of communication which could be disturbed by the use of such enzymes as well. This should be avoided or at least kept on a low level. The optimal control allows for an optimal use of the enzymes, where and when they are most helpful.
The precise mechanism of signaling molecules degradation in bacterial biofilms due to external enzyme impact remains to be elucidated. One aspect is for example the reduced diffusivity for the autoinducers as well as the applied enzymes for the control. Further numerical simulations will focus on dependence of the strategy of external addition of enzymes on the spatial organization of the bacterial communities. The latter requires consideration of an at least 2D bacterial quorum sensing model and the introduction of nonhomogeneous spatial dynamics of bacterial populations into this model, to have it even more realistic. Nevertheless, the study provides first insights for a best-possible treatment by such enzymes and how to set up a useful mathematical framework which can be further adapted and extended.