and computer research of a nonlinear stochastic models describing the dynamics of interacting populations

Abstract

The construction of nonlinear three-dimensional models of interconnected communities number dynamics is considered, taking into account competition in populations of victims. A qualitative research of the systems is carried out, equilibrium states are found, the species number dynamics graphs are constructed. For these models, an estimate of the model parameters is given and local phase portraits are constructed. The transition to the corresponding stochastic models is made. In stochastic cases, the method of constructing self-consistent stochastic models is used. A comparative analysis of deterministic and stochastic models is carried out. Efects typical for three-dimensional models with regard to competition in prey populations are revealed. A software package for the numerical solution of diferential equations systems by modified Runge-Kutta methods is used as a software tool for researching the model. The software package allows performing numerical experiments based on the implementation of algorithms for generating trajectories of multidimensional Wiener processes and multipoint distributions and algorithms for solving stochastic diferential equations. The formulation of the optimal control problem is proposed. Computer research of the models makes it possible to obtain the results of numerical experiments on the search for trajectories and the estimation of parameters. The results obtained can find application in problems of ecological systems computer modeling, as well as in problems of synthesis, optimal control and analysis of the multidimensional stochastic models stability describing the dynamics of interacting populations. computer modeling, nonlinear model of population dynamics, optimal control, stochastization of one-step processes, symbolic computation libraries, software package htp:/ceur-ws.org CEUR Workshop Proceedings (CEUR-WS.org) IS N1613-073 Woodstock'21: Symposium on the irreproducible science, June 07-11, 2021, Woodstock, NY

-

CEUR Workshop Proceedings

1. Introduction

An important tool in solving problems of predicting the state of natural systems and managing them is mathematical modeling. To solve these problems, both traditional and new methods and approaches are used. Ecological systems with various interconnections between subsystems and a change in the structure of these interconnections in the process of functioning leads to the mathematical models construction, the analytical research of which is very dificult.

Mathematical models of the interacting communities dynamics taking into account competition and with food chains are considered in [1, 2, 3, 4, 5, 6, 7]. In [2], the stability conditions of the «predator–two preys» system are investigated. In [3], a mathematical model of a system with two competing prey and one predator is analyzed, the influence of predation on the species coexistence is described. A model of two competitors’ prey dynamics with the addition of a predator species to change the competition results is studied in [4]. In [5], the deterministic stability of the three-dimensional model «predator–two preys» limit cycles is investigated. In [6, 7], three-dimensional models of population dynamics with competition and with trophic chains are considered.

In the process of stochastic modeling for various dynamical systems a method for constructing self-consistent one-step models [8] is proposed and a software package [ 9] is developed. Some systems of population dynamics based on the construction of stochastic self-consistent models are considered in [10, 11, 12, 13, 14]. In [12, 13, 14, 15, 16], a number of control problems for the models of population dynamics are considered.

When modeling population systems, various software tools are used that provide ample opportunities for conducting computational experiments. In [17, 18], the models research is carried out using the Python language and symbolic computation libraries.

In this paper we considered several types of three-dimensional models, taking into account the competition among prey species and predation are studied. The stability of these models is investigated, the equilibrium states are calculated. The transition from deterministic models to stochastic ones is performed. A computer research is carried out to study the stability. The estimation of the model parameters is carried out, the phase portraits of the system are constructed, as well as the graphs of the population size dynamics in the deterministic and stochastic cases. The research is carried out using a software package for constructing stochastic dynamic models and searching for appropriate trajectories, as well as the Jupyter application package. The obtained efects are analyzed. The formulations of optimal control problems for the models with trophic chains are proposed. 2. The deterministic models description We consider a model described by diferential equations of the form 1̇ = 1( 1 − 1 − 2 − 1 ), 2̇ = 2( 2 − 1 − 2 − 2 ), =̇ (− + 1 1 + 2 2 − ), where 1 is the population density of the first competitor, 2 is the population density of the second competitor, is the population density of the predator, is the reproduction rate of the competitor’s population in the absence of a predator, is the specific rate of consumption of the prey population by the predator population, di is the conversion factor of the prey biomass consumed by the predator in own biomass, = 1, 2 , is the coeficient of intraspecific competition, is the coeficient of interspecific competition, is the natural mortality of the predator, is the intraspecific competition of the predator, (0) ≥ 0 , (0) ≥ 0, = 1, 2 . According to the ecological sense, the constraints on the coeficients are: > 0 , ≥ 0 , ≥ 0 , > 0, > 0, > 0, > 0 , = 1, 2 . 22 = , 12 = 21 = .

Model ( 1 ) is a modification of the model described in [ 1] with the following notation: 11 = We introduce the following notation: = 1 2 + 2 1, = 1 1 + 2 2, = ( 2 − 2) − + , 1∗ = 2∗ = ∗ = ( 2 1 − 1 2) 2 + ( 2 − 1) + ( 2 − 1) ( 2 − 2) + − ( 2 − 2) + − ( 1 2 − 2 1) 1 + ( 1 − 2) + ( 1 − 2) ( 2 − 2) + ( 1 2 + 2 1) − ( 1 1 + 2 2).

( 2 − 2) + − , ,

Equilibrium states of the ( 1 ) in general form are found. The indicated equilibrium states are as follows: 0(0, 0, 0), 1 (0, 0, − ) , 2 (0, 2 , 0) , 3 ( 1 , 0, 0) ,

, 2 + 2 2 2 − + 2 2 + 2 2 ) , 5 ( + 1 1

, 0, 1 + 1 1 1 − + 1 1 ) , 6 ( 2 − 1 1 − 2 , 2 − 2 2 − 2

, 0) , 7 ( 1∗, 2∗, ∗).

The state of equilibrium 7 is an internal state of equilibrium for which the condition of positivity is satisfied. Permanent coexistence of populations in the model ( 1 ) is established under the following conditions: 1) > 0 , ≥ 0 , ≥ 0 , > 0, > 0, > 0, > 0 , = 1, 2 ; 2) there is a unique internal equilibrium state 7 such that ≠ 0 and 1∗ > 0, 2∗ > 0, ∗ > 0; 4) at least one of the following inequalities holds 1 > 2 or 2 > 1 .

Next, we consider the model: 1̇ = 1( − 11 1 − 12 2 − ), 2̇ = 2( − 21 1 − 22 2 − ), =̇ ( + 1 + 2 − ), where for = = 1 and = = 2 are the coeficients of intraspecific competition, for not equal to are the coeficients of interspecies competition. According to the ecological sense, the constraints on the coeficients are: ≥ 0, ≠ , > 0, = , , = 1, 2 , ≥ 0 , > 0 , > 0 , > 0 , > 0 . The meaning of other parameters included in the system ( 2 ) is similar to the model ( 1 ).

Model ( 2 ) is a modification of the model described in [ 1] with the following values: 1 = 2 = , 1 = 2 = , 1 = 2 = .

Next, we introduce the following notation:

= ( 11 − 12 − 21 + 22), = ( 11 22 − 12 21), = + , ̂1∗ = ( + )( 22 − 12), ̂2∗ = ( + )( 11 − 21), ̂ ∗ = − .

Equilibrium states of the model ( 2 ) in general form are found. The indicated equilibrium states are as follows:

6 ( 11 +

, 0, ) , 5 ( − 11 + , 0) , 3 ( 11

, 0, 0) , ( 22 − 12) ( 11 − 21) ,

, 0) , 11 , ) , 7 ( ̂1∗, ̂2∗, ̂ ∗).

methodology. system in general form

The state of equilibrium 7 is an internal state of equilibrium for which the condition of positivity is satisfied. Permanent coexistence of populations in the model ( 1 ) is established under the following conditions: 1) 12 ≥ 0, 21 ≥ 0, 11 > 0, 22 > 0, ≥ 0 , > 0 , > 0 , > 0 , > 0 ; 2) there is a unique internal equilibrium state 7 such that ≠ 0 and ̂1∗ > 0, ̂2∗ > 0, ̂ ∗ > 0; 4) at least one of the following inequalities holds 22 > 12 or 11 > 21.

In [1], a theoretical research of the generalized model «predator–two preys» stability is carried out. However, the numerical research of this model in order to identify the conditions of oscillating modes causes a number of dificulties. We carry out a computer research of the models ( 1 ) and ( 2 ), which are special cases of the model [1]. Subsequently, the transition to the stochastic case is made based on the method of self-consistent stochastic models.

3. Transition to stochastic models

We carry out the transition to stochastization of the models ( 1 ) and ( 2 ) using the method of constructing self-consistent stochastic models. This method is based on a combinatorial We write down the scheme of interaction between elements for the «predator–two preys»

19–32 ( 3 ) the main ones. description:

In the interaction scheme ( 3 ), the first line corresponds to the natural reproduction of prey in the absence of other factors, the second line symbolizes intraspecific (at = ) interspecific (at ≠ ) competition, and the third describes the predator-prey relationship. The fourth line is responsible for intraspecific competition among predators, and the fith describes their natural mortality.

The method of constructing self-consistent stochastic models assumes, in the course of mathematical transformations, a transition from the interaction scheme to obtaining the coeficients of the Fokker–Planck equation. This transition is carried out using the upgraded software package described in [9]. This software package is implemented in the Python programming language using the NumPy and SyPy libraries.

The software package consists of the following modules: IS_to_SDE.py and stochastic.py. The IS_to_SDE.py module is designed to obtain the coeficients of the Fokker–Planck equation from the interaction scheme. The stochastic.py module is a module for obtaining solutions for the stochastic model.

The algorithm of the software package is shown in Diagram 1.

The IS_to_SDE.py module takes as input the matrices M and N of the system states before and after interaction, vectors K_plus and K_minus interaction coeficients) and a vector X (the system state vector). As a result, we obtain a symbolic representation of the Fokker–Planck equation coeficients. Using the SymPy library allows you to get the code of these coeficients in TeX, which makes them easier to read.

The IS_to_SDE.py module consists of several functions. Hereafter there is a description of The S_plus function for obtaining the forward interaction coeficients has the following def S_plus(X, K_plus, M): """ interaction coefficient [s^{-}_{1}(x),...,s^{-}_{s}(x)]""" res = [] for i in range(len(K_plus)): return res

Prod_s = [Prod_(x, int(n)) for (x, n) in zip(X, M[i, :])] res.append(K_plus[i]* sp.prod(Prod_s))

The derivation function drift_vector for obtaining the drift vector A in the Fokker–Planck equation is described as follows:

In addition, we use the following description of the diffusion_matrix function to obtain the difusion matrix B in the Fokker–Planck equation: def diffusion_matrix(X, K_plus, K_minus, N, M): """ diffusion matrix""" res = sp.zeros(rows=len(X), cols=len(X)) R = M.T - N.T R = sp.Matrix(R) for i in range(len(K_plus)): res += R[:, i] * R[:, i].T * S(X, K_plus, K_minus, N, M)[i] return res

In fig. 2. the result of the functions for obtaining the Fokker–Planck equation coeficients for the model ( 1 ) is presented.

Further, the obtained coeficients are transferred to the stochastic.py module for the numerical solution of the generated stochastic diferential equations and drawing the graphs of these solutions.

The stochastic.py module can be used to study and numerically solve systems of ordinary diferential equations and their corresponding stochastic diferential equations based on the Runge–Kutta method and its modifications. A detailed description of this module is performed in [8, 9].

4. Results of computer experiments

For the models ( 1 ) and ( 2 ), we carry out a series of computational experiments using the above-described software package designed to study and numerically solve systems of ordinary diferential equations and the corresponding stochastic diferential equations. Computational experiments are carried out for both the deterministic case and the stochastic case.

For the model ( 1 ), calculations are carried out at = . We consider the following sets of parameters: ( 1, 2, ) = (2, 1.6, 1.2), ( 1, 2, , , 1, 2, 1, 2, ) = (0.2, 0.4, 0.8, 0, 0.2, 0.4, 0.3, 0.6, 1). With this set of parameters, the approximate equilibrium states are obtained.

Further, for the model ( 1 ), we consider the following sets of parameters: ( 1, 2, ) = (0.5, 0.4, 0.3), ( 1, 2, , , 1, 2, 1, 2, ) = (0.2, 0.4, 0.8, 0.2, 0.2, 0.4, 1.3, 2.4, 0.1). With this set of parameters, we obtained the approximate equilibrium states.

Figures 4 and 5 show the dynamics of population density for given sets of initial conditions. The dashed line indicates the fluctuations of species number of the deterministic model ( 1 ), the solid line indicates the stochastic one. Taking into account the results presented in Fig. 4, we note that the trajectories have an oscillating character with conservation of amplitudes. For the corresponding stochastic model, damping of oscillations takes place with an approximation to the stationary mode. Taking into account the results presented in Fig. 5, we note the proximity of the trajectories of the deterministic and stochastic models. In the case under consideration, stochastization does not afect the of the system behavior, which is characterized by damping of oscillations.

For the model ( 2 ), the following sets of parameters are considered: ( 1, 2, ) = ( 2, 1, 6 ), (, , 11, 12, 21, 22, , , ) = (8, 2.5/3, 8, 4, 4.125, 1, 1, 1, 0) . With this set of parameters, the approximate equilibrium states are obtained.

Fig. 6 shows the population density dynamics for two sets of initial conditions indicated above. For stochastic and deterministic cases, the trajectories are located close to each other. As in the deterministic case, the mean values graphs of various realizations of the stochastic diferential equation solutions reach a stationary mode.

Computational experiments show that the character of the systems ( 1 ) and ( 2 ) stability is significantly influenced by the coeficients of intraspecific and interspecific competition. Oscillations are formed if = = 0 and if 1 = 2, 1 = 2 at = . At = ≠ 0 , the oscillations have a damping character. At 1 = 2, 1 = 2 at ≠ one of the prey populations dies out. Next, we proceed to the consideration of controlled models and formulate the optimal control problem.

5. The problem of optimal control

Let us formulate an optimal control problem for a three-dimensional model of the interconnected communities number dynamics, taking into account competition in populations of preys. For a three-dimensional model ( 1 ), the controlled model is given by a system of diferential equations where = () are control functions. The input parameters is explained in section 2.

Constraints for model ( 4 ) are specified in the form 1(0) = 10, 2(0) = 20, 3(0) = 30, 1( ) = 11, 2( ) = 21, 3( ) = 31, ∈ [0, ], 0 ≤ 1 ≤ 11, 0 ≤ 2 ≤ 21, 0 ≤ 3 ≤ 31, ∈ [0, ].

With regard to problem ( 4 )–( 6 ), we consider the functional to be minimized in the form () =

3 ∫ ∑ (). 0 =1 ( 4 ) ( 5 ) ( 6 ) ( 7 )

Control quality criterion ( 4 ) corresponds to minimizing losses from regulation of the population, and in this case, the positive coeficients are denoted by in ( 7 ).

For the model ( 4 ), the optimal control problem can be formulated as follows: find the minimum of functional ( 7 ) under conditions ( 5 ), ( 6 ) taking into account ≥ 0.

We also generalize model ( 2 ) for the controlled case and formulate the corresponding optimal control problem. At the same time, we consider the criterion of control quality, which also consists in minimizing losses from regulation of the population. necessary:

It is possible to construct the control laws 1, 2, 3 by diferent methods. For example, it is possible to use PID controllers [19] or controllers using sliding mode [20].

For a population model with competition and migration flows in [ 13], the authors use a polynomial control of the form

In this case, the model parameters are the coeficients 1 , 2 , ..., of polynomial functions. Methods of global parametric optimization [21, 22] are usually used to calculate the parameters.

We propose to use control methods based on machine learning and regulators using artificial intelligence. In particular, it is possible to use machine learning in combination with controllers based on fuzzy logic or artificial neural networks [ 23, 24]. The generalized algorithm for constructing the optimal trajectory of the model ( 4 ) based on machine learning is shown in

Thus, the optimal control problem is to find = () those that satisfy, firstly, the phase constraints of problem ( 5 ), and secondly, the optimality criterion ( 7 ). To solve the problem, it is

a) to construct the loss function,

b) to build a parametric control model,

c) to use the global parametric optimization algorithm for search the coeficients of the control model with the minimum loss function.

We propose the construction of stochastic controlled models taking into account the «predator–two preys» interaction. To study such models, it is advisable to consider the control laws 1, 2, 3 using the algorithm for constructing the optimal trajectory of the model ( 4 ). For population models with model migration, a number of computer experiments are carried out in [12, 14]. In [14], the «predator–prey» model is studied taking into account migration lfows. In [ 12], multidimensional models of competing populations with migration, without trophic chains, are studied. In [13, 14], the indicated approach to the construction of stochastic controlled models is efective. A similar approach can be used for the models ( 1 ) and ( 2 ).

6. Conclusion

Computer research of two competing individuals interaction of prey with a predator population nonlinear models makes it possible to study the proposed models stability. The obtained results of the models research at diferent variables sets and initial conditions make it possible to estimate the influence of the predator species on the result of the prey competition. It is established that the presence of trophic chains has a positive character on the result of competition, and this, in turn, contributes to the coexistence of species. By the aid of developed software package, graphs of population dynamics are constructed. For the systems ( 1 ) and ( 2 ), the transition to the corresponding stochastic diferential equations is carried out. The introduction of stochastics makes it possible to take into account the probabilistic nature of the reproduction processes and species death, as well as random fluctuations occurring in the environment in time and leading to random fluctuations of the model parameters. We formulated an optimal control problem,

References

[1] V. Hutson , G. T. Vickers , A criterion for permanent coexistence of species, with an application to a two-prey one-predator system , Mathematical Biosciences 63 ( 1983 ) 253 - 269 . doi: 10 .1016/ 0025 - 5564 ( 82 ) 90042 - 6 .

[2] K. Fujii , Complexity-stability relationship of two-prey-one-predator species system model: Local and global stability , Journal of Theoretical Biology 69 ( 1977 ) 613 - 623 . doi: 10 .1016/ 0022 - 5193 ( 77 ) 90370 - 8 .

[3] S. B. Hsu , Predator-mediated coexistence and extinction , Mathematical Biosciences 54 ( 1981 ) 231 - 248 . doi: 10 .1016/ 0025 - 5564 ( 81 ) 90088 - 2 .

[4] M. van Baalen , V. Křivan , P. C. J. van Rijn , M. W. Sabelis , Alternative food, switching predators, and the persistence of predator‐prey systems , The American Naturalist 157 ( 2001 ) 512 - 524 . doi: 10 .1086/319933.

[5] L. V. Bashkirceva , I. A . andKarpenko, L. B. Ryashko , Stochastic sensitivity of limit cycles for «predator-two preys» model , Izvestiya VUZ. Applied Nonlinear Dynamics 18 ( 2010 ) 42 - 64 . doi: 10 .18500/ 0869 - 6632- 2010- 18- 6- 42 - 64.

[6] A. V. Epifanov , V. G. Tsibulin , On the dynamics of symmetric systems , Computer research and modeling 9 ( 2017 ) 799 - 813 .

[7] E. A. Aponina , Y. M. Aponin , A. D. Bazykin , Analysis of complex dynamic behavior in the predator-two prey model, Problems of ecological monitoring and modeling of ecosystems 5 ( 1982 ) 163 - 180 .

[8] M. N. Gevorkyan , T. R. Velieva , A. V. Korolkova , D. S. Kulyabov , L. A. Sevastyanov , Stochastic runge-kutta software package for stochastic diferential equations , in: W. Zamojski, J. Mazurkiewicz , J. Sugier , T. Walkowiak , J. Kacprzyk (Eds.), Dependability Engineering and Complex Systems , Springer International Publishing, Cham, 2016 , pp. 169 - 179 . doi: 10 .1007/978- 3- 319 - 39639- 2_ 15 .

[9] M. Gevorkyan , A. Demidova , T. Velieva , A. Korolkova , D. Kulyabov , L. Sevastyanov, Implementing a method for stochastization of one-step processes in a computer algebra system , Programming and Computer Software 44 ( 2018 ) 86 - 93 . doi: 10 .1134/S0361768818020044.

[10] A. V. Demidova , Equations of population dynamics in the form of stochastic diferential equations , RUDN Journal of Mathematics, Information Sciences and Physics ( 2013 ) 67 - 76 .

[11] O. V. Druzhinina , O. N. Masina , Methods of stability research and controllability of fuzzy and stochastic dynamical systems , Dorodnitsyn Computing Center of RAS, Moscow, 2009 .

[12] A. V. Demidova , O. V. Druzhinina , O. N. Masina , E. D. Tarova , Computer research of nonlinear stochastic models with migration flows , in: Proceedings of the Selected Papers of the 9th International Conference “Information and Telecommunication Technologies and Mathematical Modeling of High-Tech Systems” (ITTMM-2019). CEUR Workshop Proceedings , volume 2407 , 2019 , pp. 26 - 37 .

[13] A. V. Demidova , O. V. Druzhinina , O. N. Masina , A. A. Petrov , Computer research of the controlled models with migration flows , in: Proceedings of the Selected Papers of the 10th International Conference “Information and Telecommunication Technologies and Mathematical Modeling of High-Tech Systems” (ITTMM-2020). CEUR Workshop Proceedings , volume 2639 , 2020 , pp. 117 - 129 .

[14] A. V. Demidova , O. V. Druzhinina , M. Jacimovic , O. N. Masina , N. Mijajlovic , N. Olenev , A. A. Petrov , The generalized algorithms of global parametric optimization and stochastization for dynamical models of interconnected populations , in: N. Olenev , Y. Evtushenko , M. Khachay , V. Malkova (Eds.), Optimization and Applications , volume 12422 , Springer International Publishing, Cham, 2020 , pp. 40 - 54 . doi: 10 .1007/978- 3- 030 - 62867- 3 _ 4 .

[15] O. Kuzenkov , G. Kuzenkova, Optimal control of auto-reproduction systems, Bulletin of the Russian Academy of Sciences, Theory and control systems ( 2012 ) 26 - 37 .

[16] I. Danilova , The problem of optimal behavior of the population in the area taking into account intraspecific competition and migration, in: Proceedings of the XIII All-Russian Meeting on Control Problems of the VSPU-2019 , ICS of RAS, Moscow, 2019 , pp. 1476 - 1479 .

[17] C. Fuhrer , J. Solem , O. Verdier , Scientific Computing with Python 3 , Packt Publishing , 2016 .

[18] R. Lamy , Instant SymPy Starter, Packt Publishing, 2013 .

[19] E. Nikulin , Fundamentals of the theory of automatic control. frequency methods of analysis and synthesis of systems, in: Textbook for universities, BHV-Petersburg, Saint- Petersburg , 2004 .

[20] S. Ulyanov , N. Nefedov , A fuzzy controller with a sliding mode based on soft calculations, Software products and systems , Software & Systems 4 ( 2014 ) 174 - 177 . doi: 10 .15827/ 0236 - 235X . 108 .

[21] R. Storn , K. Price , Diferential evolution - a simple and eficient heuristic for global optimization over continuous spaces , Journal of Global Optimization 11 ( 1997 ) 341 - 359 . doi: 10 .1023/A: 1008202821328 .

[22] J. Lampinen , A constraint handling approach for the diferential evolution algorithm , in: Proceedings of the 2002 Congress on Evolutionary Computation. CEC'02 (Cat. No.02TH8600) , volume 2 , 2002 , pp. 1468 - 1473 . doi: 10 .1109/CEC. 2002 . 1004459 .

[23] A. Uskov , A. Kuzmin , Intelligent control technologies, Artificial neural networks and fuzzy logic, Hot line- Telecom , Moscow, 2004 .

[24] D. Rutkovskaya , M. Pilinsky , L. Rutkovsky , Neural networks, genetic algorithms and fuzzy systems, Hot line- Telecom , Moscow, 2006 .