On Targeted Control over Trajectories of Dynamical Systems Arising in Models of Complex Networks

: The question of targeted control over trajectories of systems of differential equations encountered in the theory of genetic and neural networks is considered. Examples are given of transferring trajectories corresponding to network states from the basin of attraction of one attractor to the basin of attraction of the target attractor. This article considers a system of ordinary differential equations that arises in the theory of gene networks. Each trajectory describes the current and future states of the network. The question of the possibility of reorienting a given trajectory from the initial state to the assigned attractor is considered. This implies an only partial control of the network. The difﬁculty lies in the selection of parameters, the change of which leads to the goal. Similar problems arise when modeling the response of the body’s gene networks to serious diseases (e.g., leukemia). Solving such problems is the ﬁrst step in the process of applying mathematical methods in medicine and pharmacology.


Introduction
Let us start with the following citation: "Complex systems are networks made of a number of components that interact with each other, typically in a nonlinear fashion.Complex systems may arise and evolve through self-organization, such that they are neither completely regular nor completely random, permitting the development of emergent behavior at macroscopic scales.Complex systems science is a rapidly growing scientific research area that fills the huge gap between the two traditional views that consider systems made of either completely independent or completely coupled components. . .These properties can be found in many real-world systems, e.g., gene regulatory networks within a cell, physiological systems of an organism, brains and other neural systems" [1].
We proceed with considering gene regulatory networks (GRNs in short).Living cells in an organism form complicated systems that can be studied using mathematical methods as well.The aim of these studies is to understand the complexity of these systems and the structure of their interrelations.Every element of such systems can influence others activating or inhibiting them.The gene regulatory system (GRN) is defined as a network of genes and their activating-inhibiting connections.Different mathematical models were used to analyze networks [2].Models using differential equations are especially effective since they treat networks as dynamical objects and involve the concept of an attractor.Differential equations allow for describing oscillatory behavior, stationary solutions, and cyclical patterns.Nonlinear ordinary differential equations are widespread mathematical tools for studying the regulatory interactions between genes.The time-dependent variables x(t) represent the concentration of gene products mRNAs or proteins.These variables are positive-valued.
It was noticed by biologists that cells of living organisms are adaptable to changes in an environment even if these changes are very rapid.It was proposed to use attractor selection as the principal mechanism of adaptation to unknown changes in biological systems [3].
The main idea of attractor selection is that the system is driven by deterministic and stochastic components.Attractors are a part of the solution space.Conditions of such system are controlled by a simple feedback.When conditions of a system are suitable (close to one of the attractors), it is driven almost only by deterministic behavior, and stochastic influence is very limited.When conditions of the systems are poor, the system is driven mostly by stochastic behavior.In this case, the system randomly fluctuates in search for a new attractor.When it is found, deterministic behavior again dominates over stochastic [4].
On the other hand, the system can be controlled by changing the adjustable parameters (if any).Then, stochastic behavior can be neglected (this is our assumption) and only the deterministic model can be studied.If we use the attractor selection mechanism for network resource management, at first we should define a regulatory matrix W = {w ij }, which shows relationships between node pairs, that is, how each node pair affects each other including itself.As it was proposed by some authors ( [5], for instance), three types of influence exist, namely activation, inhibition, and no relation, corresponding to positive, negative values of w ij in the interval [−1, 1], or zero.We do not restrict the range of values for the entries w ij .
Some authors consider GRNs in the conditions of serious disease [6][7][8].The mathematical model consists of a system of ordinary differential equations, which possess some remarkable properties.This system depends on multiple parameters, some of which are adjustable.The properties of this system imply the existence of attractors.These attractors can also be multiple.The above-mentioned authors associate the disease with special states of GRNs.It is claimed that the trajectory, which reflects the current system state, will tend to a "wrong" attractor.This can be improved by reasonably selected control means.Mathematically, these means are imagined as the changing of the adjustable parameters so that the trajectory changes its direction and goes to an attractor corresponding to a normal state.
In this article, we consider models of GRN, consisting of ordinary differential equations.We elaborate the scheme of control and managing trajectories in a GRN network.The aim is to redirect a trajectory from an initial point to a targeted attractor.Examples for two-dimensional systems are provided.
The problem of treating complex networks, and modeling them using mathematical means and notions, is very important due to the existence of networks in nature and technology.In our reference list, we have collected some sources, which are useful for first addressing the problem.In [9], the main objectives for the study on the border of biology and mathematics were discussed.As one of the main problems, the understanding of "the structure and the dynamics of the complex intercellular web of interactions that contribute to the structure and function of a living cell" was manifested.In [10], configurations of networks are discussed, focusing on links and nodes overlapping and considering mostly physical networks.In [11], the notion of "sensors" is introduced.It is noted that, in most cases, not all parameters can be treated explicitly, and principles of management of networks should be invented making use of only a group of available parameters.A notion of an observable system is invented.An important problem of estimating the internal state of a system from experimentally available data is discussed.In [12], the controllability of complex networks is discussed.It starts with the declaration that "the ultimate proof of our understanding of natural. . .systems is reflected in our ability to control them".Let us mention several remarks.It is noted that "a framework to control complex self-organized systems is lacking".It is known that genetic networks in a model are driven by systems with sparse regulatory matrices W. The authors of [12] have considered this point.They conclude that "sparse inhomogeneous networks, which emerge in many real complex systems, are the most difficult to control".In many sources, the controllability of a system is explained as "a dynamical system is controllable if, with a suitable choice of inputs, it can be driven from any initial state to any desired final state in a finite time".In our article, the goal is moderate.We wish to indicate means that will help to redirect a trajectory from an initial position in the phase space to a desired attractor, which is usually a stable equilibrium.The book [13] in Part II provides a great amount of information on the topic of Control of Nonlinear Systems.Several models of Control Design are proposed.Available methods of Nonlinear Control Design are discussed, as well as Robust and Adaptive controls, mathematical tools for control, and much more.
When considering genetic systems, the knowledge of the structure of phase space and the influence of parameters on the phase space structure increase the effectiveness of mathematical modeling significantly.We suggest that the geometrical approach, based on the study of isoclines and their locations, is rather natural and may lead to deeply penetrating into the essence of the problems' results.We try to illustrate this point by our treatment of a two-dimensional case.The results for higher dimensional systems need more facts and examples.In the reference list, however, one can find articles, concerning genetic systems of orders three, four, six [8,14,15], and even general ones, for arbitrary n [16,17].
We also consider systems of ordinary differential equations, which appear in the neurodynamics theory (we call them ANN systems, from Artificial Neural Networks).These systems naturally research in parallel to GRN systems, since both types of systems have many similarities.There are, however, essential differences, such as lacking parameters and the broader region of action in ANN differential equations, which have a similar structure and exhibit similar behavior.The parameters and their meaning are different, however.Therefore, the treatment of ANN model needs special attention.It is to be mentioned that neural networks are typically not associated directly with differential equations, but with difference equations or maps.
We focus on problems of control and management of GRN and ANN systems.
In order to become familiar with the topic, the resources [9-13,16-22] are useful.

GRN System
The dynamical system of the form is used to model genetic regulatory networks [2,5] and telecommunications networks [4] as well.This system first appeared in [23].The function f (z) is a sigmoid function, that is, monotonically increasing from 0 to 1 as z changes from −∞ to +∞, having only one point of inflexion, like the function 1 1 + e −µ z , v g is a parameter that controls deterministic behavior and η is stochastic term.We neglect the stochastic term in (1), so η = 0. Neglecting the stochastic term and assuming v g = 1, θ i = θ for all i, we can write the dynamical system in extended form where w ij are entries of the regulatory matrix W.
The equilibrium states can be detected from the system ( The current state of the system is described by the vector x(t).Attractors of systems of the form (2) were studied in [16,17].
General properties of systems (2) were studied and the results can be found in the related literature [2,5,21].
Of the most importance to our analysis of these systems are two facts.
The unity cube Q n = {0 < x i < 1, i = 1, . . ., n} is an invariant domain for systems of the form (2). The vector field on the border of Q n is directed inside.This can be understood by the inspection of the vector field on the faces of Q n .The difference f i (. ..) − x i is either positive or negative, depending on the choice of a face of Q n .
The second remarkable fact about systems of the form ( 2) is that their nullclines, defined by (3), can intersect only within Q n .They do, and at least one equilibrium exists.The total number of equilibria depends on the dimensionality and parameters of a system, but it is finite.
Multiple critical points of different nature can occur, depending on the choice of parameters µ and θ and elements of the regulatory matrix W.Even for n = 2 and simple W, the number of isolated critical points can be up to nine.

Attractors
We will denote the attractors of the system (4) A i .Each attractor has its basin of attraction, denoted B i .Each B i is a subset of the phase space (x 1 , . . ., x n ).If the current system state X(t) is in B i , then the system state vector X(t) will tend to B i .Attractors different of the equilibrium points are also available.Periodic attractors can be constructed easily.Examples of periodic attractors for 2D, 3D, and 4D GRN systems can be found in [14], as well as the discussion and related references.Periodic solutions in GRN systems with steep sigmoid functions were studied in [24].Chaotic attractors can appear in GRN systems, but examples are rare.
Attractors can also be distinguished by the property to be "undesired" and "normal".In real substances, this may correspond to the disease state of an organism and, respectively, to the healthy state [7].Therefore, the problem of driving the system from the undesired state (that is, in the basin of attraction of some equilibrium) to a normal state arises.This is the problem of the controllability type that is generally difficult to solve.We propose the schemes of how to drive the system to a normal state.We will also show how these schemes work in a particular situation.This particularity is due to the specific regulatory matrix W, which corresponds to the case of the cross-activation network.The respective regulatory matrix is of the form The system (4) then takes the form

Influence of Parameters on the Structure of the Phase Plane of 2D GRN Systems
The general system of ordinary differential equations which is often used to model genetic networks, in case of the two dimensions (two-element network) is It is a quasi-linear system, where the nonlinearity is represented by the logistic function f (z) = 1/(1 + e −µz ).Parameters µ and v are positive.Let us describe the influence of parameters on the phase plane, and especially on the mutual location of isoclines.Isoclines are curves, where x 1 or x 2 are constant.Especially useful are nullclines, which are given by the relations (9) Critical points (alternatively, equilibria) are solutions of the system of two equations (9).Let us list the properties of the system (8).Some of these properties are evident.Proofs of other properties are scattered over the related literature, mentioned above.

1.
There is an invariant set with the following properties: 1a The vector field defined by the system ( 8) is directed inward on the border of Q 2 ; 1b The nullclines (9) can intersect only in Q 2 ; 1c The nullclines intersect at least once.The total number of intersections is finite.For the 2D case, the maximal number of critical points is nine.For this, both nullclines have to be Z-shaped;

2.
By changing θ i , the nullclines can be shifted; by changing θ 1 , the first nullcline can be shifted in the vertical direction; by changing θ 2 , the second nullcline is moved horizontally, preserving shape; 3.
By changing µ i , the shapes of the nullclines can be changed; for sufficiently large values of µ i , the three segments of a sigmoid curve, representing a nullcline, become almost straight.In this case, the system ( 8) is almost piece-wise linear; for the study of the case of piece-wise linear system consult [24]; 4.
By changing v i , the parallelepiped Q 2 can be made stretched or compressed; for More on the classification of systems by properties of the regulatory matrices can be found in [19,25].
Remark 1.The system (8) with the matrix where a > 0, b > 0 can have one, two, or three critical points.This depends on other parameters.The case µ 1 = mu 2 = µ, θ 1 = θ 2 = theta was studied, and the region Ω was defined in the (µ, θ)-plane, which decomposes the plane with respect to the number of critical points.
Remark 2. The system (8) with the matrix where a > 0, k > 0 can have one stable critical point; then, (under k increasing) a stable periodic trajectory, and then multiple critical points, of which some are attractive.
Remark 3. The conditions for the system (8) to have a single critical point were obtained in [26].If this point is non-attractive (a saddle, or a repelling one), then the system has a limit cycle (through Andronov-Hopf bifurcation).

Inhibition Case in 2D GRN Systems
Consider the two-dimensional (2D abbreviated) system of ODE of the form (2) Look at Figures 1-3.Calculations are performed and pictures are created using Wolfram Mathematica tools, see Appendix A. Let the regularity matrix in (13) be This corresponds to the inhibition case.The nullclines (red and black) intersect three times.The green circle in Figure 1 corresponds to the current state of the 2D system.Due to the vector field, the current state is in the basin of attraction of the lower critical point, which is a stable node.The phase plane for system (13) with the matrix ( 14), µ 1 = µ 2 = 8, θ 1 = 0.01, θ 2 = −0.9.

Controllability by Changing θ
The goal is to redirect the current trajectory, emanating from the green spot, to the upper right critical point, which is conventionally the "normal" one.This can be achieved by manipulating the adjustable parameter θ.Change θ 1 from its current value −0.5 to the value 0.1.This corresponds to the shift of the first nullcline (black one) down.As the result, only one, the upper left, critical point remains, and their type is not changed.It is a stable node.The effect of this action is seen in Figure 2. The flow of the vector field will lead the green spot to the left upper, now unique, critical point, which is identified as "normal" attractor.The goal is achieved, and the system will go to the right state.

Controllability by Changing Both θ
It is clear that changing both parameters θ in (13) will lead to the shifting of both nullclines.The second nullcline (red one) can move in a horizontal direction.The change of the parameters θ 1 and θ 2 from their current values to the values 0.01 and −0.9, respectively, will lead to the configuration depicted in Figure 3.The selected trajectory will go to the desired attractive critical point at the upper-left corner.Proposition 1.In the case of inhibition (the regulatory matrix is ( 14)), any of the side critical points can be made a unique global attractor by appropriate choices of the parameters θ.

Systems of the form
arise in the theory of artificial neural networks ( [27], Chapter 6).The hyperbolic tangent function tanh z is a sigmoid function, but its range of values is (−1, 1).The invariant domain for the system ( 15) is the open cube G n = {−1 < x i < 1, i = 1, 2, . . ., n}.The nullclines are defined by the equations The cross-points of the nullclines are critical points (equilibria).At least one critical point exists in G n for the n-dimensional system (15).There is much similarity between GRN systems and ANN systems.
Consider the two-dimensional version of ( 15) The nullclines of the system (17) are defined by the equations Our goal is to establish the control over the ANN system.This system has less parameters than the GRN system.The only possibility for the system (17) to be controlled changing the parameters is to change the entries of the matrix We will show that this is possible.Look at Figures 4 and 5. Let the regularity matrix in (13) be  The nullclines (red and black) intersect three times, as seen in Figure 4. Suppose that the trajectory, corresponding to the current system state, tends to the upper-right critical point.It is a stable node.

Controllability by Changing an Element of A
The goal is to redirect the current trajectory to the lower-left critical point, which is also attractive.For this, we can only change some entries of the matrix A. Let the new matrix be Figure 5 shows a new configuration of nullclines.There is a unique critical point of the type stable node.The trajectory will go to the desired attractor.The goal is achieved.
Consider a more complicated case.Let the coefficients of the system (17) be the entries of the matrix The nullclines and the vector field are depicted in Figure 6.There are nine critical points, of which four are attractive.The green spots stand for the initial states of the system (17).The trajectories starting from these points will go to the stable critical points at the upper-left and lower-right locations.Let them be conventionally "undesired" ones.The goal is to redirect them to the attractive critical points at the upper-right and lower-left positions.This can be achieved by manipulating the elements of the matrix A. Let the element "6" in (22) be changed to the value "4".The new matrix is The new configuration of nullclines is depicted in Figure 7.The green spots are (approximately) on the border of the basins of attraction of the desired critical points.When released, the trajectories will go to the new attractors at the upper-right and/or lower-left locations, depending on where they are exactly, in the basin of attraction of the upper attractor, or in the basin of attraction of the opposite one.The goal is achieved.Proposition 2. The trajectories of the system (17) can be redirected from a given attractive critical point to another one by changing the elements of the matrix A.

Conclusions
Control over GRN systems and ANN systems is possible if by control we mean changing the properties of a system in the desired direction.In particular, it can be implemented by manipulating the nullclines.This is easier for GRN systems since they have more parameters.The most promising and geometrically understandable is changing the θ parameters.In ANN systems, the nullclines can be manipulated by the elements of the matrix A. Knowledge of the basins of attraction is a prerequisite for the implementation of control.The bistable GRN system of differential equations modeling the activation case or inhibition case can be driven from one attractor to another using several techniques.First, elements of the regulatory matrix W can be changed appropriately.Second, the parameter θ can control this process.
The following quote outlines possible further research in this direction."Because of the conceptual similarities between engineering and biological regulatory mechanisms, . . .these tools are now being used to analyze biochemical and genetic networks" [28] (p. 1).

Figure 5 .
Figure 5.The phase plane for the system (17) with the matrix (21).

Figure 6 .Figure 7 .
Figure6.The phase plane for the system(17) with the matrix(22), where the green circles denote the initial states.