Intervention in Biological Phenomena via Feedback Linearization.

The problems of modeling and intervention of biological phenomena have captured the interest of many researchers in the past few decades. The aim of the therapeutic intervention strategies is to move an undesirable state of a diseased network towards a more desirable one. Such an objective can be achieved by the application of drugs to act on some genes/metabolites that experience the undesirable behavior. For the purpose of design and analysis of intervention strategies, mathematical models that can capture the complex dynamics of the biological systems are needed. S-systems, which offer a good compromise between accuracy and mathematical flexibility, are a promising framework for modeling the dynamical behavior of biological phenomena. Due to the complex nonlinear dynamics of the biological phenomena represented by S-systems, nonlinear intervention schemes are needed to cope with the complexity of the nonlinear S-system models. Here, we present an intervention technique based on feedback linearization for biological phenomena modeled by S-systems. This technique is based on perfect knowledge of the S-system model. The proposed intervention technique is applied to the glycolytic-glycogenolytic pathway, and simulation results presented demonstrate the effectiveness of the proposed technique.

1. Introduction

Biological systems are complex processes with nonlinear dynamics. S-systems are proposed in [1, 2] as a canonical nonlinear model to capture the dynamical behavior of a large class of biological phenomena [3, 4]. They are characterized by a good tradeoff between accuracy and mathematical flexibility [5]. In this modeling approach, nonlinear systems are approximated by products of power-law functions which are derived from multivariate linearization in logarithmic coordinates. It has been shown that this type of representation is a valid description of biological processes in a variety of settings. S-systems have been proposed in the literature to mathematically capture the behavior of genetic regulatory networks [6-13]. Moreover, the problem of estimating the S-system model parameters, the rate coefficients and the kinetic orders, has been addressed by several researchers [12, 14–16]. In [17], the authors studied the controllability of S-systems based on feedback linearization approach.

Recently, the authors in [18] developed two different intervention strategies, namely, indirect and direct, for biological phenomena modeled by S-systems. The goal of these intervention strategies is to transfer the target variables from an initial steady-state level to a desired final one by manipulating the control variables. The complexity of the nonlinear biological models led researchers to focus on nonlinear control approaches, such as sliding mode control that was introduced in [19].

A basic problem in control theory is how to use feedback in order to modify the original internal dynamics of nonlinear systems to achieve some prescribed behavior [20]. In particular, feedback linearization can be used for the purpose of imposing, on the associated closed-loop system, a desired behavior of some prescribed autonomous linear system. When the system to be controlled is linear time-invariant system, this is known as the problem of pole placement, while in the more general case of nonlinear systems, this is known as the problem of feedback linearization [21, 22]. Significant advances have been made in the theory of nonlinear state feedback control, such as feedback linearization and input-output decoupling techniques [21, 22]. The state feedback linearization technique has been widely utilized in many applications. For example, the authors in [23] have used feedback linearization in cancer therapy, where full knowledge of the state and parameter vectors is assumed to transform a multiinput multioutput nonlinear system into a linear and controllable one using nonlinear state feedback. Then, linear control techniques can be applied for the resulting system [22, 24].

Hence, in this paper we consider the problem designing a nonlinear intervention strategy based on feedback linearization for biological phenomena modeled by S-systems. In this proposed algorithm, the control variables are designed such that an integral action is added to the system. The main advantage of the integral action is in improving the steady state performance of the closed-loop system. As a case study, the proposed intervention strategy is applied to a glycolytic-glycogenolytic pathway model. The glycolytic-glycogenolytic pathway model is selected as it plays an important role in cellular energy generation when the level of glucose in the blood is low (fasting state) and glycogen has to be broken down to provide the substrate to run glycolysis. By controlling the glycogenolytic reaction, one can exert control over whether glycolysis will run or not under low-glucose conditions.

This paper is organized as follows. In Section 2, the S-system model is presented and the control problem is formulated. In Section 3, some mathematical preliminaries as well as the feedback linearizable control scheme are presented. In Section 4, the glycolytic-glycogenolytic pathway model is considered as a case study. Finally, concluding remarks and possible future research directions are outlined in Section 5.

2. S-System Presentation and Problem Formulation

Consider the following S-system model [25]: where α > 0 and β > 0 are rate coefficients and θ and μ are kinetic orders and there exist N + m variables (genes/metabolites) where the first N variables are dependent and the remaining m variables are independent variables. Assume that p out of the N dependent variables are target (or output) variables (i.e., genes/metabolites that need to be regulated to some desired final values), where these output variables are defined as and i ∈ 𝒴 ⊂ {1,…, N}, where 𝒴 is the set of indices corresponding to the dependent variables that are selected as output variables. The steady-state analysis of the S-system model [1, 18] shows that when the number of dependent variables with prespecified desired values is equal to the number of independent variables (which means that we have enough degrees of freedom), the above S-system model equations will have a unique steady-state solution under the nonsingularity assumption. Hence, in order to control the expressions/concentrations of the target variables, we consider an integral control approach where the following r equations are added to the above S-system: where i ∈ 𝒰 ⊂ {N + 1,…, N + m}, where 𝒰 is the set of indices corresponding to the independent variables that are used as control variables. This means that r out of the m independent variables will be used as control variables, and the overall system will have p inputs and p outputs. It should be noted that the formulation above can be easily extended to deal with systems having more inputs than outputs. Let us denote by 𝒳 = {1,…, N} + 𝒰, where 𝒳 corresponds to the indices of all variables except the independent variables that are not used as control variables. Here, it is assumed that the values (expressions/concentrations) of the independent variables that are not used as control variables are known constants (i.e., x = δ , i ∈ {N + 1,…, N + m} − 𝒰, where δ are known constants) [6].

Figure 1 shows the S-system (1) augmented by the integral control. The S-system with integral control (1)–(3) can be written in the form where x = [x ] ∈ ℝ, i ∈ 𝒳, u = [u 1,…, u ] ∈ ℝ, y = [y 1,…, y ] ∈ ℝ and which can be expressed as where g (x) = [01,02,…,0,1,0,…,0], for i = 1,…, p.

Figure 1

S-system with integral control architecture.

Problem Formulation

Suppose that the outputs of the S-system (1) are initially at the steady-state condition y 0, j = 1,…, p. Let us denote by y , j = 1,…, p, the desired final steady state values of the output (target) variables. Then, the main goal of the feedback linearizable controller is to determinate the control inputs u , j = 1,…, p, that can guide the target variables from the initial steady-state condition to the final one [18].

3. Feedback Linearizable Intervention

Here, we show how feedback linearization can be utilized to design a nonlinear intervention strategy to control biological phenomena modeled by S-systems. Feedback linearization can be used to obtain a linear relationship between the output vector y and a new input vector v, by making a right choice of the linearizing law. Once the equivalent model becomes linear, we may design a dynamic control law-based classical linear control theory. Before starting the development of this control technique, it is important to introduce the following mathematical preliminaries [20-22].

3.1. Mathematical Preliminaries

Let the vector function f : ℜ → ℜ be a vector field in ℜ . The vector function f(x) is called a smooth vector function if it has continuous partial derivatives of any required order [26]. Given a scalar function h(x) and a vector field f(x), we define a new scalar function L h, called the Lie derivative of h with respect to f, as follows.

Definition 1 (see [26])

Let h : ℜ → ℜ be a smooth scalar function, and f : ℜ → ℜ be a smooth vector field on ℜ , then the Lie derivative of h with respect to f is a scalar function defined by L h = ∇hf.

Thus, the Lie derivative L h is simply the directional derivative of h along the direction of the vector f. Repeated Lie derivatives can be defined recursively as follows: Similarly, if g is another vector field, then the scalar function L L h(x) can be described as

Definition 2 (see [26])

Let f and g be two vector fields on ℜ . The Lie bracket of f and g is a third vector field defined by where the Lie bracket [f, g] is commonly written as ad g (where ad stands for “adjoint”).

Repeated Lie brackets can then be defined recursively by ad (0) g = g,…, ad ( g = [f, ad ( g].

3.2. Feedback Linearizable Controller

Consider the S-system model (6). Differentiating the jth output y of this system with respect to time, we get for j = 1,2, 3,…p. Note in (7) that if each of the L h (x) = 0, then the inputs do not appear in the equation. Define γ to be the smallest integer such that at least one of the inputs appears in y (, that is with at least one of the L (L ( h ) ≠ 0, for some x. Let the p × p matrix D(x) be defined as Based on the above definitions, the relative degree for multiinput multioutput (MIMO) systems is defined next.

Definition 3 (see [27])

The system (6)-(7) is said to have vector relative degree γ 1, γ 2,…, γ at x 0 if L L ( h (x) ≡ 0, 0 ≤ k ≤ γ − 2, for i = 1,…, p and the matrix D(x 0) is nonsingular.

If a system has well-defined vector relative degree, then (12) can be expressed as where Since D(x 0) is nonsingular, it follows that D(x) ∈ ℜ is bounded away from nonsingularity for x ∈ U, a neighborhood U of x 0. Then, the state feedback control law yields the linear closed-loop system The block diagram of the linearized system is shown in Figure 2.

Figure 2

Diagram block of the linearizable system.

Feedback linearization transforms the system into a linear system where linear control approaches can be applied. Here, v represents the new input vector of the linearized system.

In the case the system has vector relative degree, where γ 1 + ⋯+γ = n, the nonlinear system can be converted into a controllable linear system, where the feedback control law is defined in (16) and the coordinate transformation is ξ(x) = [L ( h (x)], 0 ≤ j ≤ γ − 1, 0 ≤ i ≤ p. Let the following distributions be defined as [27] for i = 1,…, n − 1, then we have the following result.

Proposition 4 (see [27])

Suppose that the matrix g(x 0) has rank p. Then, there exist p functions λ 1,…, λ , such that the system has vector relative degree (γ 1,…, γ ) with γ 1 + γ 2 + ⋯+γ = n if

for each 0 ≤

the dimension

for each 0 ≤

The proof of this proposition can be found in [27].

The new control vector v = [v 1,…,v ] is designed based on the desired closed-loop response, which can be written as for j = 1,…, p, where {y , y (1),…, y (, y (} denotes the desired reference trajectories for the outputs. The proportional gains are chosen such that the following polynomial is a Hurwitz polynomial [28]: The block diagram of the closed-loop system in the feedback linearizable form is shown in Figure 3.

Figure 3

Closed loop of the linearizable system.

4. Case Study

In this section, we demonstrate the efficacy of the feedback linearizable intervention approach described in this paper by applying it to a well-studied biological pathway model representing the glycolytic-glycogenolytic pathway shown in Figure 4 [17, 29]. Glycolysis is the process of breaking up a six-carbon glucose molecule into two molecules of a three-carbon compound, and glycogenolysis is the process by which the stored glycogen in the body is broken up to meet the needs for glucose. In glycogenolysis, the phosphorylase enzyme acts on the polysaccharide glycogen to reduce its length by one glucose unit. The glucose unit is released as a glucose-1 phosphate. The glycolytic-glycogenolytic pathway can be mathematically represented by the following S-system model:

Figure 4

Glycolytic-glycogenolytic pathway [29].

In this case, N = 3, m = 7 and the parameter are defined as α 1 = 0.077884314, θ 14 = 0.66, θ 16 = 1, β 1 = 1.06270825, μ 11 = 1.53, μ 12 = −0.59, μ 17 = 1, α 2 = 0.585012402, θ 21 = 0.95, θ 22 = −0.41, θ 25 = 0.32, θ 27 = 0.62, θ 210 = 0.38, β 2 = α 3 = 0.0007934561, μ 22 = θ 32 = 3.97, μ 23 = θ 33 = −3.06, μ 28 = θ 38 = 1, β 3 = 1.05880847, μ 33 = 0.3, and μ 39 = 1. Here, the model variables are defined as follows: x 1 is glucose-1-P, x 2 is glucose-6-P, x 3 is fructose-6-P, x 4 is inorganic phosphate ion, x 5 is glucose, x 6 is phosphorylase a, x 7 is phosphoglucomutase, x 8 is phosphoglucose isomerase, x 9 is phosphofructokinase, and x 10 is glucokinase.

For this model, the metabolites x 4 through x 10 are defined as independent variables, which are the variables that are not affected by other variables, and the metabolites x 1 through x 3 are defined as the dependent variables, which are the primary variables of interest that we wish to control. Here, we choose the independent variables x 4, x 5, and x 8 as manipulated or control variables, as shown in Figure 4, as they can affect the production of the dependent variables x 1, x 2, and x 3. Also, we choose to keep the independent variables x 6, x 7, x 9, and x 10 fixed ignoring their effect on the controlled variables, and assuming that the controller only uses the independent variables x 4, x 5, and x 8 to control the dependent variables x 1, x 2, and x 3. The independent variables have the following values x 4 = 10, x 5 = 5, x 6 = 3, x 7 = 40, x 8 = 136, x 9 = 2.86, and x 10 = 4. Here, we try to control x 1, x 2, and x 3 by manipulating x 4, x 5, and x 8, so we have and all other x ′s for i = 6,7, 9, and 10 are kept fixed. The initial values of the outputs y 1, y 2, and y 3 are selected as 0.067, 0.465, and 0.150, respectively, and the desired reference outputs are selected as y = 0.2, y = 0.5, and y = 0.4.

Hence, the overall system can be expressed in the form of (6), where where a 1 = α 1 x 6 , a 2 = α 2 x 7 x 10 , a 3 = α 3, b 1 = β 1 x 7 , b 2 = β 2, and b 3 = β 3 x 9 .

Based on the S-system model describing the glycolytic-glycogenolytic pathway, it can be verified that the outputs need to be differentiated twice with respect to time so that the input variables (u 1, u 2, or u 3) appear in the expressions of differentiated outputs, as follows: where Hence, in this case the system has vector relative degree γ = [γ 1,γ 2,γ 3] = [2,2,2], and hence we have γ 1 + γ 2 + γ 3 = 6.

The matrix form of the system of differential equations presented in (25) can be written in the form of (14), where The matrix D(x) is invertible if the following condition is satisfied: Based on (25), it can be seen that the control variables u 1 and u 2 appear only in the expressions of y 1 (2) and y 2 (2), respectively. However, u 3 appears in the expressions of y 2 (2) and y 3 (2). Hence, u 1 and u 3 need to be used to control y 1 and y 3, respectively, and both u 2 and u 3 are needed to control y 2.

Hence, the control laws based on (16) can be expressed as Substituting the expressions of the control variables (29) in (25), we obtain the following decoupled linear system: The new control variables v , for j = 1,2, 3, need to be designed so that the target variables y track some desired reference trajectories, y .

Using (20), the new control variables v , for j = 1,2, 3, are found to be The new control components, v 1, v 2, and v 3, are defined in (31), where the parameters are selected as k 1 = 1, k 11 = 5, k 2 = 10−3, k 21 = 20, k 3 = 3, and k 31 = 5.

Figures 5 and 6 show the output response and the control input signals when the feedback linearizable controller is applied. It is clear from Figure 5 that the system outputs converge to their desired values. Another simulation study is implemented for a different reference trajectory, where the value of the reference signal increases linearly before saturating at the desired final value. The closed-loop output response in this case is shown in Figure 7 and the control signals are shown in Figure 8. It is clear from Figure 7 that the feedback linearizable controller is driving the target variables to track the desired reference trajectories.

Figure 5

Closed-loop outputs for constant reference signals.

Figure 6

Control signals for constant reference signals.

Figure 7

Output response for closed-loop tracking.

Figure 8

Control signals for closed-loop tracking.

To study the robustness properties of the feedback linearizable controller, similar simulation studies have been conducted when the parameters μ 22 and β 2 are varied within 10% of their nominal values. It has been found that the closed-loop system is stable only for parameter variations within 1% and with unacceptable performance. This agrees with our earlier assumption that full system knowledge is needed for proper operation of the feedback linearizable controller.

5. Conclusion

In this paper, feedback linearizable control has been applied for intervention of biological phenomena modeled in the S-system framework. As a case study, the glycogenolytic-glycolytic pathway model has been used to demonstrate the efficacy of feedback linearization in controlling biological phenomena modeled by S-system. One main drawback of this approach is that it assumes full knowledge of the biological system model. Usually, the S-system model does not perfectly represent the actual dynamics of the biological phenomena. Hence, one future research direction is to develop an adaptive intervention strategy that is capable of controlling the biological system even in the presence of model uncertainties. Another future research direction is to develop intervention techniques that take into account additional constraints due to the nature of the drug injection process. Definitely, incorporating such knowledge from medical practitioners would require imposing constraints on the magnitude, duration, and possibly the rate of change of the injected drug into the design of intervention technique.