Probiotic therapy involves using live microorganisms, primarily bacteria and yeasts, to improve or restore the balance of beneficial bacteria in the body, particularly in the gut. These microorganisms, when administered in adequate amounts, can offer health benefits to the host. Probiotic therapy is used for various conditions, including diarrhea, irritable bowel syndrome (IBS), and even to support immune function. The dynamics of probiotic therapy is extremely nonlinear. Bifurcation analysis is a powerful mathematical tool used to deal with the nonlinear dynamics of any process. Several factors must be considered, and multiple objectives must be met simultaneously.  Bifurcation analysis and multi-objective nonlinear model predictive control (MNLMPC) calculations are performed on two dynamic models of probiotic therapy. The MATLAB program Matcont was used to perform the bifurcation analysis. The MNLMPC calculations were performed using the optimization language Pyomo in conjunction with the state-of-the-art global optimization solvers IPOPT and BARON. The bifurcation analysis revealed the existence of branch points in both models. The branch points (which cause multiple steady-state solutions from a singular point) are very beneficial because they enable the Multiobjective nonlinear model predictive control calculations to converge to the Utopia point (the best possible solution) in the models. It is proved (with computational validation) that the branch points were caused because of the existence of two distinct separable functions in one of the equations in each dynamic model. A theorem was developed to demonstrate this fact for any dynamic model.

Mattar et al (2001) [1] studied the effect of probiotics on enterocyte bacterial translocation in vitro Dani et al (2002) [2] studied the use of Probiotics feeding in the prevention of urinary tract infections. Millar et al (2003) [3] investigated the use of probiotics for preterm infants. Bin-Nun et al (2005) [4] studied the use of oral probiotics to prevent necrotizing enterocolitis. Land et al (2005) [5] showed that Lactobacillus sepsis was associated with probiotic therapy. Szajewska et al (2006) [6] investigated the efficacy of probiotics in gastrointestinal diseases in children. Hammerman and Kaplan (2006) [7] discussed the connection between probiotics and neonatal intestinal infection. Barclay et al (2007) [8] reviewed the use of probiotics for necrotizing enterocolitis. AlFaleh et al (2008) [9] studied the use of probiotics for the prevention of necrotizing enterocolitis in preterm infants. Lin et al (2008) [10] showed that oral probiotics prevent necrotizing enterocolitis in very low birth weight preterm infants. Claud and Walker (2008) [11] studied bacterial colonization, probiotics, and necrotizing enterocolitis. Arciero et al (2010) [12] developed a mathematical model to Analyze the Role of Probiotics and Inflammation in Necrotizing Enterocolitis. Zhang et al (2015) [13] investigated the impacts of gut bacteria on human health and diseases. Ahmed and Jawad (2023) [14] performed a bifurcation analysis of the role of good and bad bacteria in the decomposing toxins in the intestine with the impact of antibiotic and probiotics supplement. This work aims to perform bifurcation analysis and multi-objective nonlinear control (MNLMPC) studies in two models involving probiotics, which are discussed in Arciero et al (2010) [12] (model 1), and Ahmed and Jawad (2023) [14] (model 2). The paper is organized as follows. First, the model equations are presented, followed by a discussion of the numerical techniques involving bifurcation analysis and multi-objective nonlinear model predictive control (MNLMPC). The results are then presented, followed by the discussion and conclusions.

Model Description

Probiotic Model 1Arciero et al (2010) [12]

The model equations are 

Here, bl represents the pathogenic bacteria in the intestinal lumen, bpbl represents the Probiotic bacteria in the intestinal lumen, the permeability of the intestinal wall to bacteria, b is the pathogenic bacteria in the blood/tissue, bpb represents the probiotic bacteria in the blood/tissue, and mv represents the activated inflammatory cells.  

The parameter values are 

r2 was used as the bifurcation parameter and the control value.

Probiotic model 2Ahmed and Jawad (2023) [14]

The dynamic model equations are 

(b1, b2, c, a) represent the good bacteria, the bad bacteria, the non-decomposing toxins in the large intestine, and the concentration of dissolved antibiotics. The base parameter values are.

r1 was used as the bifurcation parameter and control value.

Bifurcation analysis 

The MATLAB software MATCONT is used to perform the bifurcation calculations. Bifurcation analysis deals with multiple steady-states and limit cycles. Multiple steady states occur because of the existence of branch and limit points. Hopf bifurcation points cause limit cycles. A commonly used MATLAB program that locates limit points, branch points, and Hopf bifurcation points is MATCONT (Dhooge Govearts, and Kuznetsov, 2003[15]; Dhooge Govearts, Kuznetsov, Mestrom and Riet, 2004[16]). This program detects Limit points (LP), branch points (BP), and Hopf bifurcation points(H) for an ODE system 

 Let the bifurcation parameter be   Since the gradient is orthogonal to the tangent vector, 

The tangent plane at any point must satisfy 

Where A is 

where is the Jacobian matrix. For both limit and branch points, the matrix  must be singular. The n+1 th component of the tangent vector = 0 for a limit point (LP)and for a branch point (BP) the matrix must be singular. At a Hopf bifurcation point, 

@ Indicates the bialternate product while is the n-square identity matrix. Hopf bifurcations cause limit cycles and should be eliminated because limit cycles make optimization and control tasks very difficult.  More details can be found in Kuznetsov (1998 [17]; 2009 [18]) and Govaerts [2000] [19]. 

Hopf bifurcations cause unwanted oscillatory behavior and limit cycles. The tanh activation function (where a control value u is replaced by) is commonly used in neural nets (Dubey et al 2022[20]; Kamalov et al, 2021 [21] and Szandała, 2020 [22) and optimal control problems (Sridhar 202[23]) to eliminate spikes in the optimal control profile. Hopf bifurcation points cause oscillatory behavior. Oscillations are similar to spikes, and the results in Sridhar (2024) [24] demonstrate that the tanh factor also eliminates the Hopf bifurcation by preventing the occurrence of oscillations. Sridhar (2024) [24] explained with several examples how the activation factor involving the tanh function successfully eliminates the limit cycle causing Hopf bifurcation points. This was because the tanh function increases the time period of the oscillatory behavior, which occurs in the form of a limit cycle caused by Hopf bifurcations.

Multi-objective Nonlinear Model Predictive Control (MNLMPC) 

Flores Tlacuahuaz et al (2012) [25] developed a multiobjective nonlinear model predictive control (MNLMPC) method that is rigorous and does not involve weighting functions or additional constraints. This procedure is used for performing the MNLMPC calculations Here (j=1, 2.n) represents the variables that need to be minimized/maximized simultaneously for a problem involving a set of ODE 

                                                     (7)

  being the final time value, and n the total number of objective variables and. u the control parameter. This MNLMPC procedure first solves the single objective optimal control problem independently optimizing each of the variables individually. The minimization/maximization of will lead to the values  .  Then the optimization problem that will be solved is 

This will provide the values of u at various times. The first obtained control value of u is implemented and the rest are discarded. This procedure is repeated until the implemented and the first obtained control values are the same or if the Utopia point where ( for all j) is obtained. 

Pyomo (Hart et al, 2017) [26] is used for these calculations. Here, the differential equations are converted to a Nonlinear Program (NLP) using the orthogonal collocation method  The NLP is solved using IPOPT (Wächter And Biegler, 2006) [27]and confirmed as a global solution with BARON (Tawarmalani, M. and N. V. Sahinidis 2005) [28].

The steps of the algorithm are as follows 

1. Optimize and obtain at various time intervals ti. The subscript i is the index for each time step. 
2. Minimize and get the control values for various times.
3. Implement the first obtained control values 
4. Repeat steps 1 to 3 until there is an insignificant difference between the implemented and the first obtained value of the control variables or if the Utopia point is achieved. The Utopia point is when for all j. 

Sridhar (2024) [29] proved that the MNLMPC calculations to converge to the Utopia solution when the bifurcation analysis revealed the presence of limit and branch points. This was done by imposing the singularity condition on the co-state equation (Upreti, 2013) [30]. If the minimization of lead to the value and the minimization of lead to the value  The MNLPMC calculations will minimize the function .  The multi-objective optimal control problem is

Differentiating the objective function results in

The Utopia point requires that both and are zero. Hence

the optimal control co-state equation (Upreti; 2013) [30] is 

is the Lagrangian multiplier. is the final time.  The first term in this equation is 0 and hence 

At a limit or a branch point, for the set of ODE  is singular. Hence there are two different vectors-values for where and . In between there is a vector where . This, coupled with the boundary condition will lead to This makes the problem an unconstrained optimization problem, and the only solution is the Utopia solution. 

Probiotic model 1

When r2 was used as the bifurcation parameter a branch point was located at 

(bl, bp, bl, eps, b, bpb, mv, r2) values of (13.2359, 0, 0.1860, 0.2626, 0.1287, 0.2987, 0.1976)) This is shown in Fig. 1. For the MNLMPC calculations, were minimized individually and led to values of 20.5914 and 0.21552. r2 was the control parameter. The multiobjective optimal control problem will involve the minimization of subject to the equations governing Model 1. This led to a value of zero (the Utopia solution). The MNLMPC control value (r2) was 00.95777. Figs 2 and 3. show the various MNLMPC profiles. Fig. 4 shows the control profile of r2. This profile exhibited noise, which was remedied by using the Savitzky-Golay filter to produce the smooth version of r2 (r2sg).

Probiotic model 2

When r1 was used as the bifurcation parameter, a branch point was located at 

(b1, b2, c, a, r1) values of (0; 0; 4.000000 5.084746; 0.532881). This is shown in Fig. 5.

For the MNLMPC calculations, were minimized individually and led to values of 0, 8 and 5.0847. r1 was the control parameter.  The multiobjective optimal control problem will involve the minimization of subject to the equations governing Model 1. This led to a value of zero (the Utopia solution).   The MNLMPC control value (1) was 0.2185649. Figs 6-9 and 3. show the various MNLMPC profiles. 

                                                                             Figure 1: Bifurcation analysis Probiotic model 1 (indicating branch point)

                                                                                               Figure 2: MNLMPC Probiotic model 1 (bl, bpbl, bpb)

                                                                                      Figure 3: MNLMPC Probiotic model 1 (mv, eps,b)

                             Figure 4: MNLMPC Probiotic model 1 (r2, r2sg) (r2sg is the smooth version of r2 obtained by using the Savitzky Golay Filter)

                                                                          Figure 5: Bifurcation analysis Probiotic model 2 (indicating branch point)

                                                                                                      Figure 6: MNLMPC Probiotic model 1 (b2)

                                                                                                   Figure 7: MNLMPC Probiotic model 1 (a)

                                                                                                          Figure 8: MNLMPC Probiotic model 1 (c)

                                                                                                          Figure 9: MNLMPC Probiotic model 1 (r1)

Theorem

If one of the functions in a dynamic system is separable into two distinct functions, a branch point singularity will occur in the system. 

Proof

Consider a system of equations 

  . Defining the matrix A as 

is the bifurcation parameter. The matrix A can be written in a compact form as 

The tangent at any point x; () must satisfy 

The matrix  must be singular at both limit and branch points. The n+1 th component of the tangent vector at a limit point (LP) and for a branch point (BP) the matrix must be singular. 

Let any of the functions fi are separable into 2 functions as 

At steady-state and this will imply that either or or both and must be 0.  This implies that two branches and  will meet at a point where both and are 0. 

At this point, the matrix B will be singular as a row in this matrix would be 

The singularity in B implies that there exists a branch point.

In the probiotic model 1, a branch point was located at 

(bl, bpbl, eps, b, bpb, mv, r2) values of (13.2359, 0, 0.1860, 0.2626, 0.1287, 0.2987, 0.1976))

(Here, the two distinct functions can be obtained from the second ODE in probiotic model 1

The two distinct functions are 

and

The values of satisfy both the equations and computationally validate the theorem.

In the probiotic model 2, a branch point was located at the values (b1, b2, c, a, r1) of (0, 0, 4.000000, 5.084746, 0.532881). 

(Here, the two distinct functions can be obtained from the first ODE in probiotic model 2,                                                                                                      

The two distinct functions are 

and                                           

Setting 

satisfies both the equations and validates the theorem.

Additionally, the MNLMPC calculations in both models converge to the Utopia solution justifying the analysis of Sridhar (2024) [29]. 

Bifurcation analysis and multiobjective nonlinear control (MNLMPC) studies in two probiotic therapy models. The bifurcation analysis revealed the existence and branch points in both models. The branch points (which cause multiple steady-state solutions from a singular point) are very beneficial because they enable the Multiobjective nonlinear model predictive control calculations to converge to the Utopia point (the best possible solution) in the models. It is proved (with computational validation) that the branch points were caused because of the existence of two distinct separable functions in one of the equations in each dynamic model. A theorem was developed to demonstrate this fact for any dynamic model.  A combination of bifurcation analysis and Multiobjective Nonlinear Model Predictive Control (MNLMPC) for dynamic models involving probiotic therapy is the main contribution of this paper. 

All data used is presented in the paper

The author, Dr. Lakshmi N Sridhar, has no conflict of interest.

Dr. Sridhar thanks Dr. Carlos Ramirez and Dr. Suleiman for encouraging him to write single-author papers