=Paper= {{Paper |id=Vol-2407/paper-03-166 |storemode=property |title= Computer research of nonlinear stochastic models with migration flows |pdfUrl=https://ceur-ws.org/Vol-2407/paper-03-166.pdf |volume=Vol-2407 |authors=Anastasia V. Demidova,Olga V. Druzhinina,Olga N. Masina,Ekaterina D. Tarova |dblpUrl=https://dblp.org/rec/conf/ittmm/DemidovaDMT19 }} == Computer research of nonlinear stochastic models with migration flows == https://ceur-ws.org/Vol-2407/paper-03-166.pdf
26


UDC 004.4
        Computer research of nonlinear stochastic models with
                          migration flows
                     Anastasia V. Demidova* , Olga V. Druzhinina† ,
                        Olga N. Masina‡ , Ekaterina D. Tarova‡
                      *
                       Department of Applied Probability and Informatics
                            Peoples’ Friendship University of Russia
                Miklukho-Maklaya str. 6, Moscow, 117198, Russian Federation
                 †
                   Federal Research Center “Computer Science and Control”
                                of Russian Academy of Sciences
               Vavilov str. 44, building 2, Moscow, 119333, Russian Federation
                                ‡
                                  Bunin Yelets State University
                  Communards str. 28, Yelets, 399770, Russian Federation
     Email: demidova_av@rudn.university, ovdruzh@mail.ru, olga121@inbox.ru, katerina.tarova@yandex.ru

   The problems of the design methods extension and computer research of nondeterministic
finite-dimensional population models describing migration flows are studied. The significant
difficulties arise in the construction of high dimension dynamic models in the course of
analytical research. Computer research allows not only to obtain the results of numerical
experiments to search for trajectories and estimate the parameters of deterministic models,
but also to reveal the effects caused by stochasticization. The model parameters are estimated
and local phase portraits are constructed for the initial four-dimensional migration-population
model. The transition from the vector ordinary differential equation to the corresponding
stochastic differential equation is performed. The structure of the stochastic model is described
on the basis of applying the method of constructing self-consistent stochastic models. As a tool
for the study of population-migration models, a software package is used to solve numerically
the differential equations systems using modified Runge–Kutta methods. The software package
allows numerical experiments based on the implementation of algorithms for generating
trajectories of multidimensional Wiener processes and multipoint distributions and algorithms
for solving stochastic differential equations. The comparative analysis of the computer research
results obtained for stochastic models is carried out. The properties of migration-population
systems in deterministic and stochastic cases are characterized. The comparison of the results
obtained for the three-dimensional and four-dimensional cases is carried out. The effects
inherent in models with migration flows are revealed.The obtained results can be applied to
the problems of modeling and forecasting the behavior of multidimensional systems describing
the migration flows.

    Key words and phrases: computer modeling, nonlinear models of migration flows,
stochasticization of one-step processes, symbolic representation algorithms, symbol computing
libraries.




Copyright © 2019 for the individual papers by the papers’ authors. Copying permitted for private and
academic purposes. This volume is published and copyrighted by its editors.
In: K. E. Samouylov, L. A. Sevastianov, D. S. Kulyabov (eds.): Selected Papers of the IX Conference
“Information and Telecommunication Technologies and Mathematical Modeling of High-Tech Systems”,
Moscow, Russia, 19-Apr-2019, published at http://ceur-ws.org
                                      Demidova A. V. et al.                                    27


                                      1.   Introduction
    The construction and study of the models with migration are devoted to the work
of various authors (see, for example, [1–3]). In [2–6] migration flows in deterministic
population is considered. The analysis of distributed multidimensional population
models taking into account cross-migration is carried out in papers [7–10]. Migration
mechanisms can be described by both linear and non-linear functions, and lead to
different effects [2, 11]. In [3, 5, 11, 12] the issues of qualitative behavior and sustainability
of population-migration models are studied. In [13–15] a qualitative analysis of a
three-dimensional non-deterministic model with migration is performed. To study
the model, a combination of known methods of synthesis and analysis of models and
a method for constructing stochastic self-consistent models are used [16]. In [12] a
methodological support is developed for the analysis and synthesis of multidimensional
nonlinear dynamic models describing migration flows taking into account the effects of
broadband parametric and additive noise. The stability of stationary states is studied
and the effects obtained for stochastic models are interpreted. The model examples
show a comparison of the migration-population systems properties in deterministic and
stochastic cases and the effects due to stochastic broadband perturbations are revealed.
    When modeling population-migration systems, various software tools are used that
present wide possibilities for building computer models and carrying out computational
experiments. However, many software products do not contain libraries for numerical
and symbolic calculations and do not have sufficient computational complexity. In this
regard, in the study of the population-migration systems models, the application of
mathematical packages and general-purpose programming languages [17–19] is relevant.
In [20] sufficient conditions are proposed for the uniform stability of the equilibrium
states of a four-dimensional non-linear model of population dynamics. In [21] the
Maxima mathematical package is used to study the stability of this model. One of
the instrumental software tools for studying population-migration models is a software
package for the numerical solution of differential equations systems using modified Runge
– Kutta methods. The specified complex was developed in [22, 23]. The software package
allows numerical experiments based on the implementation of algorithms for generating
trajectories of multidimensional Wiener processes and multipoint distributions and
algorithms for solving stochastic differential equations.
    The several types of four-dimensional models with migration flows are studied in
this paper. A comparative analysis of the computer research results obtained for three-
dimensional and four-dimensional stochastic models with migration flows is carried
out. A comparison is made of the qualitative properties of four-dimensional models,
taking into account changes in migration rates, as well as intraspecific and interspecific
interaction coefficients. The properties of the models in the deterministic and stochastic
cases are characterized. As a tool for the study of the models, a software package in the
Python language using the NumPy and SciPy libraries is used.

                         2.   Deterministic Migration Models
   A nonlinear model is considered that is described by a system of ordinary differential
equations of the form:

                     𝑥˙1 = 𝑥1 − 𝑥21 − 𝑝13 𝑥1 𝑥3 − 𝑝14 𝑥1 𝑥4 + 𝛽𝑥2 − 𝛾𝑥1 ,
                     𝑥˙2 = 𝑥2 − 𝑥22 + 𝛾𝑥1 − 𝛽𝑥2 ,
                                                                                              (1)
                     𝑥˙3 = 𝑥3 − 𝑝31 𝑥1 𝑥3 − 𝑥23 − 𝑝34 𝑥3 𝑥4 ,
                     𝑥˙4 = 𝑥4 − 𝑝41 𝑥1 𝑥4 − 𝑝43 𝑥3 𝑥4 − 𝑥24 .

   The following notation is used in (1): 𝑥1 , 𝑥3 and 𝑥4 is the population density of
competing species in the area of the forma 𝑥1 (area 1), 𝑥2 is the population density
28                                                                                 ITTMM—2019


in the area of the form 𝑥2 (area 2), 𝑝13 , 𝑝14 , 𝑝31 , 𝑝34 , 𝑝41 , 𝑝43 > 0 are the coefficients of
species competition in the area of 1, 𝛽 > 0 and 𝛾 > 0 are the coefficients of the migration
of species between two areas, while area 2 is a refuge and 𝛽 ̸= 𝛾.
    In the model (1) the first two equations describe the dynamics of one type taking
into account migration processes. The first equation sets the dynamics in the first area,
the second one sets the dynamics in the second area. The third and fourth equations
describe the dynamics of the second and third types, interacting as competitors with
the first type in the first area. Previously, we studied the three-dimensional model (two
species, the first species migrated to another area, and in the first area it competed with
the second type).
    In the case when 𝑝13 = 𝑝31 , 𝑝14 = 𝑝41 , 𝑝34 = 𝑝43 , the model (1) takes view

                     𝑥˙1 = 𝑥1 − 𝑥21 − 𝑝13 𝑥1 𝑥3 − 𝑝14 𝑥1 𝑥4 + 𝛽𝑥2 − 𝛾𝑥1 ,
                     𝑥˙2 = 𝑥2 − 𝑥22 + 𝛾𝑥1 − 𝛽𝑥2 ,
                                                                                              (2)
                     𝑥˙3 = 𝑥3 − 𝑥23 − 𝑝13 𝑥1 𝑥3 − 𝑝34 𝑥3 𝑥4 ,
                     𝑥˙4 = 𝑥4 − 𝑥24 − 𝑝14 𝑥1 𝑥4 − 𝑝34 𝑥3 𝑥4 .

     At the same migration rates 𝛽 = 𝛾 = 𝜀 from (2) we get a model of the form

                     𝑥˙1 = 𝑥1 − 𝑥21 − 𝑝13 𝑥1 𝑥3 − 𝑝14 𝑥1 𝑥4 + 𝜀(𝑥2 − 𝑥1 ),
                     𝑥˙2 = 𝑥2 − 𝑥22 + 𝜀(𝑥2 − 𝑥1 ),
                                                                                              (3)
                     𝑥˙3 = 𝑥3 − 𝑥23 − 𝑝13 𝑥1 𝑥3 − 𝑝34 𝑥3 𝑥4 ,
                     𝑥˙4 = 𝑥4 − 𝑥24 − 𝑝14 𝑥1 𝑥4 − 𝑝34 𝑥3 𝑥4 .

   In this paper, deterministic models (1)–(3) and their stochastic generalizations are
studied.
   A three-dimensional model that describes two types of dynamics, the first species
migrating to another area, and in the first area competing with the second type, has the
form:

                            𝑦˙1 = 𝑦1 − 𝑦12 − 𝑝13 𝑦1 𝑦3 + 𝛽𝑦2 − 𝛾𝑦1 ,
                            𝑦˙2 = 𝑦2 − 𝑦22 + 𝛾𝑦1 − 𝛽𝑦2 ,
                                                                                              (4)
                            𝑦˙3 = 𝑦3 − 𝑦32 − 𝑝13 𝑦1 𝑦3 ,


    where 𝑦1 and 𝑦3 are the density of populations of competing species in the area of the
form 𝑦1 (area 1), 𝑦2 is the population density in the area type 𝑦2 (area 2), 𝑝13 > 0 is the
coefficient of competition of species in the area 1. The results of the analysis of the (4)
model and the generalized stochastic model are presented in [13–15]. The model (4) is a
generalization of the model considered in [3], in the case when the migration rates are
different.
    The numerical experiment for the model (4) is carried out with the help of the
developed software for the study and numerical solution of systems of ordinary and
stochastic differential equations using the Runge–Kutta method [22, 23]. The library is
written in Python using the NumPy and SciPy modules. The models (1)–(3) and their
stochastic generalizations are studied using this software package.
    For the numerical experiment of the deterministic model, the following sets of
parameters are chosen: the initial values 𝑦1 (0) = 0.5, 𝑦2 (0) = 1.0, 𝑦3 (0) = 0.5 for
the three-dimensional model and 𝑥1 (0) = 0.5, 𝑥2 (0) = 1.0, 𝑥3 (0) = 0.5,𝑥4 = 0.7 for a
                                         Demidova A. V. et al.                                             29


four-dimensional model, the values of the parameters 𝑝13 = 1.2, 𝑝14 = 0.5, 𝑝34 = 1.4,
𝛽 = 1.5, 𝛾 = 0.2.
   Fig. 1 presents solution trajectories for (2) and (4) models with the specified initial
values and parameter values in the time interval [0,40]. The effect of an additional
competitor on the migration flow is shown. The introduction of an additional competitor
corresponding to the phase variable 𝑥4 , leads to a decrease in both the density of the
migrating population in both areas and the density of the competitor corresponding to
the variable 𝑥3 .




                  Figure 1. Model solution trajectories (2) and (4) with
    (𝑦1 (0), 𝑦2 (0), 𝑦3 (0)) = (0.8, 0.1, 1.0), (𝑥1 (0), 𝑥2 (0), 𝑥3 (0), 𝑥4 (0)) = (0.8, 0.1, 1.0, 1.0),
                         𝑝13 = 1.2, 𝑝14 = 0.5, 𝑝34 = 1.4, 𝛽 = 1.5, 𝛾 = 0.2.

    With the values of the parameters 𝑝13 = 0.9, 𝑝14 = 0.4, 𝑝34 = 0.8, 𝑝34 = 0.6,
𝑝41 = 0.5, 𝑝43 = 0.7, 𝛽 = 0.25, 𝛾 = 0.3 in the model (1) there are the following equilibria:
𝑃1 (0, 0, 0, 0), 𝑃2 (0, 0, 0, 1), 𝑃3 (0, 0, 1, 0), 𝑃4 (0, 0, 0.69, 0.517), 𝑃5 (0.63, 0.949, 0.4957, 0),
𝑃6 (0.75, 0.98, 0.043, 0.595). Here 𝑃1 , 𝑃 5 are unstable nodes, 𝑃2 − 𝑃4 are saddles, 𝑃6 is a
stable node.
    Consider the models (1) and (2) with the same set of migration velocity values for these
models, but with different sets of intraspecific and interspecific interaction coefficients.
Solution trajectories for the model (1) with the values of the parameters 𝑝13 = 0.9,
𝑝14 = 0.4, 𝑝31 = 0.8, 𝑝34 = 0.6, 𝑝41 = 0.5, 𝑝43 = 0.7 and the trajectories of solutions for
the model (2) with 𝑝13 = 1, 𝑝14 = 0.6, 𝑝34 = 0.7 in the time interval [0, 25] is presented
in fig. 2. For models (1) and (2) the initial values are (𝑥1 (0), 𝑥2 (0), 𝑥3 (0), 𝑥4 (0)) =
(0.3, 0.7, 1.1, 0.8), values of migration speeds 𝛽 = 0.25, 𝛾 = 0.3. In fig. 2 a record of the
form 𝑥𝑖 (1), 𝑖 = 1, 2, 3, 4, indicates the trajectory corresponding to the phase variable 𝑥𝑖
for the model (1), a record of the form 𝑥𝑖 (2) denotes the trajectory corresponding to
the phase variable 𝑥𝑖 for the model (2).
    As a result of a comparative analysis of the solution trajectories of (1) and (2)
presented in fig. 2:
  1) in the model (2) the population of 𝑥4 is dying out, the population density of 𝑥1 is
       increasing, and the population density of 𝑥3 is decreasing;
  2) for the models (1) and (2) favorable conditions are created in the second area, and
       the population density of 𝑥2 increases;
  3) in the model (2) competition intensified, which led to the extinction of the 𝑥4
       population and the improvement of living conditions for the 𝑥1 and 𝑥3 populations.
    For the values of the parameters 𝑝13 = 1.2, 𝑝14 = 0.5, 𝑝34 = 1.4, 𝜀 = 0.3 in model (3)
there are the following equilibrium positions: 𝑃1 (0, 0, 0, 0), 𝑃2 (0, 0, 0, 1), 𝑃3 (0, 0, 1, 0),
𝑃4 (0, 0, 0.42, 0.42), 𝑃5 (1, 1, 0, 0), 𝑃6 (0.76, 0.94, 0, 0.62). Here 𝑃2 − 𝑃5 are saddles, 𝑃1 is
an unstable node, 𝑃6 is a stable node.
30                                                                                           ITTMM—2019




                  Figure 2. Model solution trajectories (1) and (2) with
               (𝑥1 (0), 𝑥2 (0), 𝑥3 (0), 𝑥4 (0)) = (0.3, 0.7, 1.1, 0.8), 𝛽 = 0.25, 𝛾 = 0.3.

    Consider the models (2) and (3) taking into account the coincidence of the intraspecific
and interspecific interaction for these models. Note that in the model (2) in the general
case 𝛽 ̸= 𝛾, and in the model (3) 𝛽 = 𝛾. Solution trajectories for models (2) and (3)
with the values of the parameters 𝑝13 = 1.2, 𝑝14 = 0.5, 𝑝34 = 1.4 and the initial values
(𝑥1 (0), 𝑥2 (0), 𝑥3 (0), 𝑥4 (0)) = (0.5, 1, 0.8, 1) in the time interval [0, 25] are presented
in fig. 3. For the model (2) the values of migration speeds 𝛽 = 1.5, 𝛾 = 0.2, for the
model (3) the value of migration speed 𝜀 = 0.3. In fig. 3 the record of the form 𝑥𝑖 (3)
denotes the trajectory corresponding to the phase variable 𝑥𝑖 for the model (3).




                       Figure 3. Model solution trajectories (2) and (3)
     (𝑥1 (0), 𝑥2 (0), 𝑥3 (0), 𝑥4 (0)) = (0.5, 1, 0.8, 1), 𝑝13 = 1.2, 𝑝14 = 0.5, 𝑝34 = 1.4, 𝛽 = 1.5,
                                             𝛾 = 0.2, 𝜀 = 0.3.

   As a result of a comparative analysis of the solution trajectories of the (2) and (3)
models shown in fig. 3, the following conclusions are obtained:
  1) in the model (3)the population density 𝑥4 increased insignificantly, and the pop-
     ulation density 𝑥1 decreased, while the population 𝑥3 in both models is rapidly
     dying out;
  2) as the migration rate increases, the model (3) in the second area creates more
     favorable conditions than in the model (2), so the population density 𝑥2 increases
     significantly.
   For a series of parameter sets, stationary states of the model (3) are found. It is
established that for the values of the parameters 𝑝13 = 1.2, 𝑝14 = 0.5, 𝑝34 = 1.4, 𝛽 = 1.5,
𝛾 = 0.2 in the model (2) the following stationary states exist: 𝑃1 (0, 0, 0, 0), 𝑃2 (0, 0, 0, 1),
                                       Demidova A. V. et al.                                      31


𝑃3 (0, 0, 1, 0), 𝑃4 (1.18, 0.29, 0, 0), 𝑃5 (0, 0, 0.42, 0.42), 𝑃6 (0.93, 0.25, 0, 0.53). In addition,
the system is linearized in the vicinity of stationary states and the stability on the
first approach is investigated. For the model (2) the projections of the phase portraits
(𝑥1 , 𝑥2 ) are shown in fig. 4. The projections of the phase portraits (𝑥1 , 𝑥2 ) in the vicinity
of equilibrium positions are shown in fig. 5. The projections of the phase portraits
(𝑥1 , 𝑥3 ) are shown in fig. 6.




                 a) 𝑃1                                                       b) 𝑃2




                 c) 𝑃3                                                       d) 𝑃4




                 d) 𝑃5                                                       d) 𝑃6

                         Figure 4. Phase portraits for the model (2).
32                                                                          ITTMM—2019




                a) 𝑃1                                               b) 𝑃2




                c) 𝑃3                                               d) 𝑃4




                d) 𝑃5                                               d) 𝑃6

     Figure 5. Phase portraits of the projections (𝑥1 , 𝑥2 ) for the model (1) in the
                     vicinity of the equilibrium positions 𝑃1 –𝑃6 .

                          3.   Stochastic Migration Models
   Questions of constructing and studying a stochastic model for a three-dimensional
system with regard to migration are considered in the works [13–15]. In this paper,
stochastization of (1)–(3) models is carried out using the method of constructing self-
consistent stochastic models [16]. The basis of this method is combinatorial methodology.
This method involves recording the system under study as an interaction scheme, i.e. a
symbolic record of all possible interactions between elements of the system. For this, the
                     Demidova A. V. et al.                 33




a) 𝑃1
                                                   b) 𝑃2




c) 𝑃3                                              d) 𝑃4




d) 𝑃5                                              d) 𝑃6

        Figure 6. Phase portraits for the model (3).
34                                                                           ITTMM—2019


system state operators and the system state change operator are used. Then we can get
the drift and diffusion coefficients for the Fokker–Planck equation, which allows to write
the equation itself and its equivalent stochastic differential equation in the Langevin
form.
    However, an increase in the dimension and complexity of the system under study leads
to the complexity of the analytical conclusion of the necessary coefficients of the Fokker–
Planck equation. As a solution to this problem, a software implementation is developed
for obtaining the coefficients of the Fokker–Planck equation from interaction schemes
using a symbolic computation system. This software implementation is a modification
of the stochastization method for one-step processes in the computer algebra system
described in [24]. This implementation is introduced as a module into a software package
developed earlier for the numerical study of deterministic and stochastic models [22].
    The following is the algorithm for obtaining the symbolic notation of a stochastic
differential equation (Algorithm 1).
    Algorithm 1. Getting the system of differential equations from the interaction
scheme
    Initial parameters: interaction scheme.
    Result: system of differential equations in the form of Langevin.
    Start:
   1. Getting the system state operators from the interaction scheme.
   2. Getting the change of the system state.
   3. Getting the transition intensities.
   4. Record the coefficients of the Fokker–Planck equation.
   5. Record the system of differential equations.
    end
    To implement the described algorithm, SymPy [17], computer computing system
is used, which is a powerful symbolic computation library for the Python language.
In addition, the output data obtained using the SymPy library can be transferred for
numerical calculations using the NumPy [19] library and SciPy [18].
    To obtain a stochastic model corresponding to the models (1)–(3), it is necessary to
write the interaction scheme, which is as follows:
                                  𝑎
                                  𝑖
                              𝑋𝑖 −→ 2𝑋𝑖 ,    𝑖 = 1, 4;
                                  𝑏𝑖
                        𝑋𝑖 + 𝑋𝑖 −→ 𝑋𝑖 ,     𝑖 = 1, 4;
                                  𝑐𝑖𝑗
                                                                                       (5)
                        𝑋𝑖 + 𝑋𝑗 −−→ 𝑋𝑗 ,     𝑖, 𝑗 = 1, 3, 4,   𝑖 ̸= 𝑗;
                                  𝑑
                                  𝑖
                              𝑋𝑖 −→ 𝑋𝑗 ,    𝑖, 𝑗 = 1, 2,   𝑖 ̸= 𝑗.

    The values of the coefficients for each of the models are given in the table 1.
    The first line corresponds to the natural reproduction of the species in the absence
of other factors, the second line symbolizes intraspecific competition, and the third one
is interspecific competition. The last line corresponds to the description of the species
migration process from one area to another. Then, expressions for the coefficients of the
Fokker–Planck equation are obtained for this interaction scheme using the developed
software (fig. 7).
    Further, the coefficients are transferred to the appropriate module of the software
complex for the numerical solution of the constructed stochastic differential equation.
For the numerical experiment of the obtained stochastic model, the parameters are
chosen, the same as for the numerical analysis of the deterministic model (3). Fig. 8
presents the trajectories of average values for 100 realizations for the specified initial
values and parameter values in the time interval [0, 40].
    A numerical experiment showed the closeness of trajectories for stochastic and
deterministic cases. Similarly to the deterministic case, the mean values trajectories
of different implementations of the stochastic differential equation solutions go to the
                                      Demidova A. V. et al.                                       35



                                                                                           Table 1
                       Coefficients for the interaction scheme (5)


                    Model (1)               Model (2)               Model (3)
           𝑎𝑖       𝑎𝑖 = 1, 𝑖 = 1, 4        𝑎𝑖 = 1, 𝑖 = 1, 4        𝑎𝑖 = 1, 𝑖 = 1, 4
           𝑏𝑖       𝑏𝑖 = 1, 𝑖 = 1, 4        𝑏𝑖 = 1, 𝑖 = 1, 4        𝑏𝑖 = 1, 𝑖 = 1, 4
           𝑐𝑖𝑗      𝑐13 = 𝑝13 , 𝑐14 =       𝑐13 = 𝑐31 = 𝑝13 ,       𝑐13 = 𝑐31 = 𝑝13 ,
                    𝑝14 , 𝑐31 = 𝑝31 ,       𝑐14 = 𝑐41 = 𝑝14 ,       𝑐14 = 𝑐41 = 𝑝14 ,
                    𝑐34 = 𝑝34 , 𝑐41 =       𝑐34 = 𝑐43 = 𝑝34         𝑐34 = 𝑐43 = 𝑝34
                    𝑝41 , 𝑐43 = 𝑝43 ,
           𝑑𝑖       𝑑1 = 𝛾, 𝑑2 = 𝛽          𝑑1 = 𝛾, 𝑑2 = 𝛽          𝑑1 = 𝑑2 = 𝜀




            Figure 7. The coefficients of the equation of Fokker–Planck




                       Figure 8. Model Solution Trajectories (5) with
 (𝑥1 (0), 𝑥2 (0), 𝑥3 (0), 𝑥4 (0)) = (0.8, 0.1, 1.0, 1.0), 𝑝13 = 1.2, 𝑝14 = 0.5, 𝑝34 = 1.4, 𝛽 = 1.5,
                                                 𝛾 = 0.2.

stationary mode, retain the character and values of the stationary state are weakly
distinguished.
    A series of numerical experiments is carried out using different sets of model parame-
ters. For each set, we analyzed the stationary states and constructed the trajectories of
solutions. Taking into account the characteristic cases, the peculiarities of the influence
36                                                                             ITTMM—2019


of an additional competitor on the dynamics of the migrating population are revealed
and the effects of stochastization of one-step processes are described.

                                    4.    Conclusions
    A computer study of nonlinear models with migration flows made it possible obtain
the results of numerical experiments on trajectory search and parameter estimation in the
case of high dimensionality of models, as well as to reveal effects due to stochasticization.
The estimation of the influence of an additional competitor type on the dynamics, stability
and level of migration in the multidimensional model is obtained. A comparative analysis
of computer simulation results showed that with an increase in the rate of migration
in the second area, more favorable conditions are created, which leads to a significant
increase in the number of the second population. At the same time, competition in
the first area is increasing, which leads to the extinction of one of the population. The
comparison of the results obtained for the three-dimensional and four-dimensional cases.
The software package developed in the Python language using the NumPy and SciPy
libraries is demonstrated sufficient efficiency for computer studies of multidimensional
nonlinear migration models. Numerical experiments using problem-oriented software
show the closeness of the trajectories types for stochastic and deterministic cases. The
results obtained can be applied to the problems of computer modeling of ecological,
demographic and socio-economic systems.

                                   Acknowledgments
   The publication has been prepared with the support of the “RUDN University
Program 5-100” (Demidova A.V., numerical analysis).

                                         References
1.   N. Tkachenko, J. D. Weissmann, W. P. Petersen, G. Lake, C. P. E. Zollikofer,
     S. Callegari, Individual-based modelling of population growth and diffusion in
     discrete time 2 (4) (2017) e0176101. doi:10.1371/journal.pone.0176101.
2.   H. C. Tuckwell, A study of some diffusion models of population growth, Theoretical
     Population Biology 5 (3) (1974) 345–357. doi:10.1016/0040-5809(74)90057-4.
3.   X.-a. Zhang, L. Chen, The linear and nonlinear diffusion of the competitive
     Lotka–Volterra model, Nonlinear Analysis 66 (2007) 2767–2776. doi:10.1016/j.na.
     2006.04.006.
4.   Y. Svirezhev, Nonlinear waves, dissipative structures and disasters in ecology,
     Moscow: Science, 1987.
5.   Z. Lu, Y. Takeuchi, Global asymptotic behavior in single-species discrete diffu-
     sion systems, Journal of Mathematical Biology 32 (1) (1993) 67–77. doi:10.1007/
     BF00160375.
6.   J. Cui, L. Chen, The effect of diffusion on the time varying logistic population
     growth, Computers & Mathematics with Applications 36 (1998) 1–9. doi:10.1016/
     S0898-1221(98)00124-2.
7.   L. Chen, A. Jüngel, Analysis of a multi-dimensional parabolic population model
     with strong cross-diffusion, SIAM J. Math. Anal. 36 (2004) 301–322. doi:10.1137/
     S0036141003427798.
8.   Z. Wen, S. Fu, Global solutions to a class of multi-species reaction-diffusion systems
     with cross-diffusions arising in population dynamics, J. Comput. Appl. Math. 230 (2)
     (2009) 34–43. doi:10.1016/j.cam.2008.10.064.
9.   N. Zamponi, A. Jüngel, Analysis of degenerate cross-diffusion population models
     with volume filling, Annales de l’Institut Henri Poincaré C, Analyse non linéaire
     34 (1) (2017) 1–29. doi:10.1016/j.anihpc.2015.08.003.
                                   Demidova A. V. et al.                                37


10. X. Chen, E. S. Daus, A. Jüngel, Global existence analysis of cross-diffusion popu-
    lation systems for multiple species, Archive for Rational Mechanics and Analysis
    227 (2) (2018) 715–747. doi:10.1007/s00205-017-1172-6.
11. L. J. S. Allen, Persistence and extinction in single-species reaction-diffusion models,
    Bulletin of Mathematical Biology 45 (2) (1983) 209–227. doi:10.1007/BF02462357.
12. I. N. Sinitsyn, O. V. Druzhinina, O. N. Masina, Analytical modeling and stability
    analysis of nonlinear broadband migration flow, Nonlinear world 16 (3) (2018) 3–16.
13. A. Demidova, O. Druzhinina, M. Jacimovic, O. Masina, Construction and analysis
    of nondeterministic models of population dynamics, Vol. 678, 2016, pp. 498–510.
    doi:10.1007/978-3-319-51917-3_43.
14. A. V. Demidova, O. V. Druzhinina, O. N. Masina, Design and stability analysis of
    nondeterministic multidimensional populations dynamics models, Vol. 1995, 2017,
    pp. 14–21.
15. A. V. Demidova, O. V. Druzhinina, M. Jacimovic, O. N. Masina, N. Mijajlovic,
    Synthesis and analysis of multidimensional mathematical models of population
    dynamics, in: 2018 10th International Congress on Ultra Modern Telecommu-
    nications and Control Systems and Workshops (ICUMT), 2018, pp. 361–366.
    doi:10.1109/ICUMT.2018.8631252.
16. A. V. Demidova, M. N. Gevorkyan, A. D. Egorov, D. S. Kulyabov, A. V. Ko-
    rolkova, S. L.A., Influence of stochastization on one–step models, RUDN Journal of
    Mathematics, Information Sciences and Physics (1) (2014) 71–85.
17. R. Lamy, Instant SymPy Starter, Packt Publishing, 2013.
18. T. E. Oliphant, Python for scientific computing, Computing in Science and Engg.
    9 (3) (2007) 10–20. doi:10.1109/MCSE.2007.58.
19. T. E. Oliphant, Guide to NumPy, 2nd Edition, CreateSpace Independent Publishing
    Platform, USA, 2015.
20. E. D. Tarova, Application the divergent method to the stability analysis of the
    dynamic population nonlinear model, Nonlinear world 15 (3) (2017) 59–64.
21. E. D. Tarova, Qualitative research of the population dynamics four-dimensional
    model with the help of lyapunov forst method, Nonlinear world 16 (3) (2018) 17–24.
22. E. G. Eferina, A. V. Korolkova, M. N. Gevorkyan, D. S. Kulyabov, L. A. Sevastyanov,
    One–step stochastic processes simulation software package, RUDN Journal of Math-
    ematics, Information Sciences and Physics (3) (2014) 46–59.
23. M. N. Gevorkyan, T. R. Velieva, A. V. Korolkova, D. S. Kulyabov, L. A. Sevastyanov,
    Stochastic runge–kutta software package for stochastic differential equations, in:
    W. Zamojski, J. Mazurkiewicz, J. Sugier, T. Walkowiak, J. Kacprzyk (Eds.), De-
    pendability Engineering and Complex Systems, Springer International Publishing,
    Cham, 2016, pp. 169–179. doi:10.1007/978-3-319-39639-2_15.
24. M. N. Gevorkyan, A. V. Demidova, T. R. Velieva, A. V. Korol’kova, D. S. Kulyabov,
    L. A. Sevast’yanov, Implementing a method for stochastization of one-step processes
    in a computer algebra system, Programming and Computer Software 44 (2) (2018)
    86–93. doi:10.1134/S0361768818020044.