# Dissipativity and Error Feedback Controller Design of Time-Delay Genetic Regulatory Networks

Authors: Yonggang Ma; Junmei Liu; Jiao Ai
Journal: Journal of Mathematics
 (2021)

---

## Abstract

Genetic regulatory networks (GRNs) play an important role in the development and evolution of the biological system. With the rapid development of DNA technology, further research on GRNs becomes possible. In this paper, we discuss a class of time-delay genetic regulatory networks with external inputs. Firstly, under some reasonable assumptions, using matrix measures, matrix norm inequalities, and Halanay inequalities, we give the global dissipative properties of the solution of the time-delay genetic regulation networks and estimate the parameter-dependent global attraction set. Secondly, an error feedback control system is designed for the time-delay genetic control networks. Furthermore, we prove that the estimation error of the model is asymptotically stable. Finally, two examples are used to illustrate the validity of the theoretical results.


---

## Body

## 1. Introduction

In the past few decades, research on gene regulatory networks’ modeling has attracted many biologists and mathematicians. There have been many research studies in this area, such as models based on Boolean networks [1], models based on Bayesian methods [2], and models based on differential equations [3, 4]. In particular, the gene regulation networks’ model based on ordinary differential equation technology is welcomed by researchers because of its simple mathematical technology and clear biological meaning.The concept of dissipative system originates from the research of control theory. In a sense, it is a generalization of Lyapunov’s stability theory. There are fruitful research results on stability issues in various research fields, such as neural network systems [5–10], biological systems [11, 12], and nonlinear systems [13, 14]. However, the study of stability requires the existence of a balance point. In some special cases, the balance point of the system may not exist. At this time, we consider the dissipative characteristics of the system. In the dissipative problem, we need to find a positive invariant set, and we only need to consider the dynamic characteristics of the system in the invariant set. So far, the dissipative characteristics of gene regulatory network models based on ordinary differential equation technology have not been seen in related publications.In addition, the control system has been widely used in all aspects of human society. In the ecosystem, a large number of phenomena can be observed through the setting of the controller, and different control schemes will produce different kinds of effects. In [15], the author gives some characteristics of the linear time-varying system under the control of time-varying feedback conjugate. In [16], the authors establish the stability criterion and saddle node bifurcation of gene regulatory networks that have nothing to do with coupling delay. In [17], the authors discuss a delayed fractional-order two-gene regulatory network model that can more accurately reflect the memory and genetic characteristics of the gene network.Based on the enlightenment obtained from the abovementioned literature research, we discuss the characteristics of time-delay gene regulation networks with external inputs. Compared with most previous related results, our novelty lies in the use of matrix measurement theory to carry out discussions. The matrix measurement strategy has the following two advantages. (1) Lyapunov function does not need to be constructed. (2) The result obtained by the matrix measurement is more accurate and less conservative than the result obtained by the matrix norm.The rest of the paper is arranged as follows. In Section2, we give the mathematical form of the gene regulatory network model and some basic assumptions. In Section 3, the global dissipativity of the model is discussed. In Section 4, we study the asymptotic stability of the error feedback control system. In Section 5, two examples are given to verify the validity of our conclusions.Notations: in this paper,Rn represents the n-dimensional Euclidean space. AT represents the transpose of matrix A. x∈Rn∖Ω means x∈Rn but x∉Ω. For real symmetric matrices A and B, A≽0 indicates that A is positive semidefinite. A≽B represents A−B≽0. λmaxA represents the maximum eigenvalue of matrix A.


## 2. Preliminaries

The gene regulatory networks’ model has the following classic form:(1)m˙it=−aimit+∑j=1nwijhjpjt−σit+qi,p˙it=−cipit+dimit−τit,i=1,2,…,n,where mit,pit∈R are the concentrations of the ith mRNA and protein, respectively, ai and ci are positive constants, representing the degradation rate of the ith mRNA and protein, respectively, di>0 is the translation rate from the ith mRNA to the ith protein, hjx=x/βjHj/1+x/βjHj is the regulatory functions of mRNA, where Hj is the Hill coefficient, βj is a scalar, σit and τit are transcriptional and translational delay with 0<σit≤σ,0<τit≤τ, and qi=∑j∈νiαij with νi is the set of all the transcription factor j which is a repressor of gene i, and(2)wij=αij,if factorjactivates genei,0,if there is no connection betweenjandi,−αij,if factorjrepresses genei.We rewrite model (1) in the following matrix form:(3)m˙t=−Amt+Bhpt−σt+q,p˙t=−Cpit+Dmt−τt,where mt=m1t,m2t,…,mntT,pt=p1t,p2t,…,pntT, A=diaga1,a2,…,an,B=wijn×n, C=diagc1,c2,…,cn,D=diagd1,d2,…,dn, σt=σ1t,σ2t,…,σntT,τt=τ1t,τ2t,…,τntT, hpt=h1p1t,h2p2t,…,hnpntT, and q=q1,q2,…,qnT.Next, we give some basic assumptions and lemmas.Assumption 1.
The activation functionhi⋅ meets the Lipschitz condition.Remark 1.
In most literature,hi⋅ satisfied 0≤hia−hib/a−b<li, for all a,b∈R,li>0,i∈1,2,…,n. The condition can be reduced to mi−≤hia−hib/a−b≤mi+,∀a,b∈R, where mi− and mi+ are constants. Furthermore, the condition is reduced to Assumption 1, i.e., there exists a constant li>0 such that, for any x,y∈R, the following inequalities hold: hix−hiy≤lix−y,li=maxmi−,mi+. Clearly, Assumption 1 is more general and less conservative.
Letm∗,p∗ be the equilibrium point of (3), and it satisfies(4)0=−Am∗+Bhp∗+q,0=−Cp∗+Dm∗.
Letxt=mt−m∗andyt=pt−p∗, and we obtain(5)x˙t=−Axt+Bfyt−σt,y˙t=−Cyt+Dxt−τt,where xt=x1t,…,xntT,yt=y1t,…,yntT,fyt=f1y1t,…,fnyntT, and fiyit=hiyit+pi∗−hipi∗.
In the field of gene therapy, the effect of external forces or stimuli on the biological system is very important. For example, mRNA vaccine is to transfer RNA to human cells after relevant modification in vitro for expression and produce protein antigen, which can lead the body to produce immune response to antigens, and then expand the immune ability of the body. Therefore, we consider the following time-delay GRNs with external inputs described by(6)x˙t=−Axt+Bfyt−σt+u,y˙t=−Cyt+Dxt−τt+v,where u=u1,u2,…,unandv=v1,v2,…,vn are control inputs for mRNA and protein, respectively.
The initial conditions of system (6) are(7)xs=φ1s,…,φns,ys=ψ1s,…,ψns,s∈−ρ,0,ρ=maxτ,σ.
By the definition offi⋅ and Assumption 1, it satisfies(8)fix≤lix.Definition 1.
If there exists a compact setΩ⊂Rn, for any ε>0 and any initial value φs,ψs∈Rn∖Ω,s∈−ρ,0, such that ∃T>0, when t≥T,xt,φ,yt,ψ⊆Ωε, where Ωε is the ε-neighborhood of Ω; then, the time-delay GRNs (6) are said to be a global dissipative system. In this case, Ω is called a globally attractive set. If, for any initial value φs,ψs∈Ω,s∈−ρ,0, there is xt,φ,yt,ψ⊆Ω,t≥0, then the set is said to be positive invariant. Obviously, the globally attractive set is positive invariant.Definition 2.
If there exists a compact setΩ⊂Rn, for any ε>0 and any initial value φs,ψs∈Rn∖Ω,s∈−ρ,0, there exist M1φ,M2ψ>0 and α,β>0 such that(9)infx¯∈Ωεxt,φ−x¯≤M1φexp−αt,infy¯∈Ωεyt,ψ−y¯≤M2ψexp−βt.
Then, the time-delay GRNs (6) are said to be a globally exponentially dissipative system.Definition 3.
For constant matrixA=aijn×n∈Rn×n, the matrix norms are given as(10)A1=maxj∑i=1naij,A2=λmaxATA,A∞=maxi∑j=1naij.
Forx∈Rn, the vector norms are defined as follows:(11)x1=∑i=1nxi,x2=∑i=1nxi21/2,A∞=max1≤i≤nxi.Definition 4.
Suppose thatA=aijn×n is a real matrix; then, the matrix measure of A is defined as(12)μpA=limh⟶0+En+hAp−1h,where En is a n×n identity matrix, p=1,2,∞, and ⋅p is the corresponding induced matrix norm.
It can be calculated directly by Definition4:(13)μ1A=maxjajj+∑i=1,i≠jnaij,μ2A=λmaxAT+A2,μ∞A=maxiaii+∑j=1,j≠inaij.Lemma 1 (see [18]).
Letr>0,ξ>0,η>0,andτ>0, if the continuous function Vt≥0,t∈R, and(14)D+Vt≤r−ξVt+ηV¯t,t≤t0.There exists a constantσ>0 such that −ξ+η≤−σ; then,(1)
WhenV¯t0>r/σ, the following formula is established:(15)Vt≤rσ+V¯t0−rσe−μ∗t−t0,t≥t0.(2)
WhenV¯t0≤r/σ, the following formula is established:(16)Vt≤rσ,t≥t0,where V¯t0=supt−τ≤s≤tVs, μ∗ is the unique positive root of the transcendental equation μ−ξ+ηeμτ=0, and D+Vt=limh⟶0¯Vt+h−Vt/h is the upper-right Dini derivative.Lemma 2 (see [19]).
For anyA,B∈Rn×nandp=1,2,∞, the matrix measure μp⋅ has the following properties: 1−Ap≤μpA≤Ap;2μpαA=αμpA,∀α>0;3μpA+B≤μpA+μpB.


## 3. Global Dissipativity

In this section, we give the global dissipative properties of the studied mathematical model.Theorem 1.
Iffi⋅ satisfies (8), there exist positive constant δ1, and the matrix measures μp⋅,p=1,2,∞, such that Π1+Π2≤−δ1<0, where Π1=maxμp−A,μp−C and Π2=maxlBp,Dp,l=max1≤i≤lli. Then, the time-delay GRNs (6) are a globally dissipative system, and the set Ω1=x,y∈Rn|xp+yp≤up+vp/δ1 is a positive invariant set and a globally attractive set of (6).Proof.
Let(17)Vt=xtp+ytp.
We calculate the upper-right hand Dini derivative ofVt along the solution of (6):(18)D+Vt≤limh⟶0+¯xt+hp−xtph+limh⟶0+¯yt+hp−ytph=limh⟶0+¯xt+hx˙t+ohp−xtph+limh⟶0+¯yt+hy˙t+ohp−ytph≤limh⟶0+¯xt+h−Axtp+ohp−xtph+limh⟶0+¯hBfyt−σtp+huph+limh⟶0+¯yt+h−Cytp+hDxt−τtp+hvp+ohp−ytph≤limh⟶0+¯En+h−Ap−1hxtp+limh⟶0+¯En+h−Cp−1hytp+Bpfyt−σtp+Dpxt−τtp+up+vp.
From (8), we can obtain(19)fyt−σtp≤lyt−σtp.
By (18) and (19), according to the definition of matrix measure, we deduce that(20)D+Vt≤μp−Axtp+μp−Cytp+lBpyt−σtp+Dpxt−τtp+up+vp.
According to Lemma1 and φsp+ψsp>up+vp/δ1,s∈−ρ,0, we have(21)Vt≤up+vpδ1+Me−μ∗t,where μ∗ is the unique positive root of equation μ+Π1+Π2eμρ=0andM=sup−ρ≤s≤0Vs−up+vp/δ1. Then, (6) is a globally dissipative system, and the globally attractive set is(22)Ω1=x,y∈Rn|xp+yp≤up+vpδ1.
Ifφsp+ψsp≤up+vp/δ1,s∈−ρ,0, by Lemma 2, we have Vt≤up+vp/δ1, for all t≥0. Therefore, the set Ω1 is also a positive invariant set of (6).Remark 2.
Obviously, from (22), we can easily get the following sets:(23)Ω1x=x∈Rn|xp≤up+vpδ1,Ω1y=y∈Rn|yp≤up+vpδ1.Corollary 1.
Letp=2, and other conditions are the same as Theorem 1. Then, the time-delay GRNs (6) are a globally exponentially dissipative system, and the globally attractive sets are Ω1x=x∈Rn|xp≤up+vp/δ1andΩ1y=y∈Rn|yp≤up+vp/δ1.Proof.
By the proof of Theorem1, if φs2+ψs2>up+vp/δ1,xt∈Rn\Ω1x, we choose x¯ as the intersection point of the vector ox⟶ and the spherical surface of Ω1x=x∈Rn|xp≤up+vp/δ1, where o is the origin of coordinates. We have(24)xt−x¯2=xt2−x¯2=xt2−u2+v2δ1.
Accordingly,(25)infx¯∈Ω1xxt−x¯2≤xt2−u2+v2δ1≤M1e−μ1∗t.
Similarly,(26)infy¯∈Ω1yyt−y¯2≤yt2−u2+v2δ1≤M2e−μ2∗t.
Thus, the time-delay GRNs (6) are a globally exponentially dissipative system and Ω1xandΩ1y are globally exponentially attractive sets of (6).
In Theorem1, parameter li is considered only as the maximum value information. If the information of each component li is utilized, the conservativeness of the result can be decreased.
Denote(27)Q˜=−A00−C,B˜=0B00,D˜=00D0,L=diagl1,l2,…,ln,I=uv,Fyt=diagf1y1ty1t,…,fnyntynt,F˜yt=000Fyt,En=diag1,…,1,G˜=En000.
Then, (6) can be written in the form of matrix blocks:(28)z˙t=Q˜zt+B˜F˜yt−σtzt−σt+D˜G˜zt−τt+I,where zt=xTt,yTtT.Theorem 2.
Iffi⋅ satisfies (8), there exist positive constant δ2, and the matrix measures μp⋅,p=1,2,∞, such that −δ2=μpQ˜+BLp+D˜G˜p<0. Then, the time-delay GRNs (6) are a globally dissipative system, and the sets,(29)Ω2x=x∈Rn|xp≤Ipδ2,Ω2y=y∈Rn|yp≤Ipδ2,are the positive invariant sets and the globally attractive sets of (6).Proof.
Consider the upper-right hand Dini derivative ofztp with respect to time along the trajectory of (6), and we obtain(30)D+ztp=limh⟶0+¯zt+hp−ztph=limh⟶0+¯zt+hQ˜zt+hB˜F˜yt−σtzt−σt+hD˜G˜zt−τt+hIp−ztph≤limh⟶0+¯E2n+hQ˜p−1hztp+B˜F˜yt−σtpzt−σtp+D˜G˜pzt−τtp+Ip≤μpQ˜ztp+B˜F˜yt−σtpzt−σtp+D˜G˜pzt−τtp+Ip.
By the definition of⋅p and (8), it is easy to acquire that B˜F˜yt−σtp≤BLp holds for p=1,∞. Meanwhile, it also holds for p=2. In fact, for LL−Fyt−σtFyt−σt≽0, we have BLL−Fyt−σtFyt−σtBT≽0, i.e., BLLBT≽BFyt−σtFyt−σtBT. By Weyl’s monotonicity principle, we can obtain(31)λmaxBLLBT≥λmaxBFyt−σtFyt−σtBT.
SinceλBLLBT=λBTLLB and(32)λBFyt−σtFyt−σtBT=λBTFyt−σtFyt−σtB.
We can obtain(33)λmaxBTLLB≥λmaxBTFyt−σtFyt−σtBand BFyt−σt2≤BL2.
Accordingly, we have(34)B˜F˜yt−σtp=0BFyt−σt00p=BFyt−σtp≤BLp.
Combining (30) and (34), we derive(35)D+ztp≤μpQ˜ztp+BLpzt−σtp+D˜G˜pzt−τtp+Ip≤μpQ˜ztp+BLp+D˜G˜psupt−ρ≤s<tzsp+Ip.
Letzt¯p=supt−ρ≤s<tzsp, and the above inequality can reduce to(36)D+ztp≤μpQ˜ztp+BLp+D˜G˜pzt¯p+Ip.
By inequalityμpQ˜+BLp+D˜G˜p<0 and Lemma 1, we obtain(37)ztp≤Ipδ2+M3e−μ3∗t,where μ3∗ is the unique positive root of equation μ+μpQ˜+BLp+D˜G˜peμρ=0andM3=sup−ρ≤s≤0zs−Ip/δ2. Similar to the proof of Theorem 1, we can obtain that the time-delay GRNs (6) are a globally dissipative system, and Ω2x=x∈Rn|xp≤Ip/δ2andΩ2y=y∈Rn|yp≤Ip/δ2 are positive invariant sets and globally attractive sets of (6).


## 4. State Estimation of Error Feedback Control System

In this section, we will discuss the global exponential stability of the following error dynamic system:(38)x˙t=−Axt+Bfyt−σt+u,y˙t=−Cyt+Dxt−τt+v,and the considered state estimator is of the following form:(39)x˜˙t=−Ax˜t+Bfy˜t−σt+u+U,y˜˙t=−Cy˜t+Dx˜t−τt+v+V,where UV be an error feedback control term with the form UV=K1e1K2e2, e1=x˜t−xtande2=y˜t−yt. K1 and K2 are feedback control gains, K1=diagk1,k2,…,knandK2=diagk1′,k2′,…,kn′. Thus, the error dynamic system is given as(40)e˙1t=−Ae1t+Bf′e2t−σt+U,e˙2t=−Ce2t+De1t−τt+V.In order to facilitate the discussion later, we write (40) for the matrix block form(41)e˙t=K˜et+B˜F′˜e2t−σtet−σt+D˜G˜et−τt,where(42)et=e1Tt,e2TtT,K˜=−A+K100−C+K2,F′e2t=diagf1′e21te21t,…,fn′e2nte2nt,F′˜e2t=000F′e2t.Theorem 3.
Iffi⋅ meets (8), there exist matrix K1andK2 and matrix measure μp⋅,p=1,2,∞, such that −μpK˜>BLp+D˜G˜p>0. Then, the error dynamic system (40) is globally exponential stable under controller UV=K1e1K2e2.Proof.
By inequality (30), we have(43)D+etp=limh⟶0+¯et+hp−etph=limh⟶0+¯et+hK˜et+hB˜F′˜e2t−σtet−σt+hD˜G˜et−τtp−etph≤limh⟶0+¯E2n+hK˜p−1hetp+B˜F′˜e2t−σtpet−σtp+D˜G˜pet−τtp≤μpK˜etp+B˜F′˜e2t−σtpet−σtp+D˜G˜pet−τtp≤μpK˜etp+BLp+D˜G˜psupt−ρ≤s<tesp.
Letm=−μpK˜andn=BLp+D˜G˜p. If m>n>0, that is, −μpK˜>BLp+D˜G˜p>0, we have(44)etp≤supt−ρ≤s<tespe−rt,where r=m−nerρ=−μpK˜−BLp+D˜G˜perρ. Therefore, et converges exponentially to zero with convergence rate of r, i.e., the error system is globally exponential stable.


## 5. Numerical Simulations

In this part, we give a pair of examples to illustrate the validity of our theoretical results.Example 1.
Consider time-delay GRNs (6) with(45)A=C=diag2,2,2,2,2,D=diag1,1,1,1,1,u=22.891.982.961.95,v=1.281.891.981.961.95,B=0−0.4−0.400−0.4000.40.400.40000.4−0.40000000.40,σt=0.3cost+0.2,τt=12.5sint+0.5,fjx=x2x2+1,j=1,2,…,n.
By computation, we can find thatlj−=0andlj+=0.65, and hence,(46)L−=diag0,0,0,0,0,L+=diag0.65,0.65,0.65,0.65,0.65.
In order to make the calculation simple, we chooseA,C,andD to be a diagonal matrix. At this time, p takes any value in 1,2,∞, and we have μpA=2,μpC=2,μpD=2,andp=1,2,∞. Since B is a nondiagonal matrix, it can be obtained through simple calculations: B1=0.4,B2=0.48,andB∞=0.4. Obviously, the parameters satisfy the conditions of Theorem 1, i.e., Π1+Π2=−1=−δ1<0. This means that system (6) is globally dissipative. From the vector norm, we get the following results:(47)u1=11.78,v1=9.06,u2=5.36,v2=4.0957,u∞=2.96,v∞=1.98.
Whenp=1, Ω1=x,y∈R5|x1+y1≤20.84. When p=2, Ω2=x,y∈R5|x2+y2≤9.4557. When p=∞, Ω∞=x,y∈R5|x∞+y∞≤20.84. Figure 1 shows the time-dependent trajectory of mRNA and protein concentrations when the initial values are x0=0.2,0.5,0.9,1.5,2.5T and y0=0.8,1,1.8,2.5,0.6T. The verification of Theorem 2 conditions is similar to the verification in Example 2, so we will not state it here.Figure 1
(a) The trajectory of the concentration of mRNAs with the parameters and initial values in Example1 over time (t). (b) The trajectory of the protein concentration with the parameters and initial values in Example 1 over time t.
(a)(b)Example 2.
In this example, we consider the error dynamic system (40), in which the parameters are listed as follows:(48)A=diag3,3,3,C=diag4,4,4,D=diag1,1,1,K1=diag1,1,1,K2=diag1,1,1,B=0−0.40.4−0.400.40−0.40.4,σt=0.3cost+0.2,τt=12.5sint+0.5,K˜=−A+K100−C+K2,and fjx is the same as Example 1. It is easy to obtain that(49)L−=0,L+=diag0.65,0.65,0.65,when p takes any one of 1,2,∞, and there are μpK˜=2. It can be directly obtained by simple matlab programming:(50)BL1=0.78,D˜G˜1=1,BL2=0.3199,D˜G˜2=0,BL∞=0,D˜G˜∞=1.
The above calculation results are easy to see that the conditions in Theorem3 are satisfied. That is, −μpK˜>BLp+D˜G˜p>0 holds. This shows that system (40) is globally exponentially stable. Figure 2 shows the trajectory of error system (40) with time when the initial values are e10=0.8,1,1.8T and e20=1.5,0.6,2.4T.Figure 2
(a) Trajectory of mRNA concentration error with parameters and initial values in Example2 versus time (t). (b) Trajectory of protein concentration error with parameters and initial values in Example 2 versus time t.
(a)(b)


## 6. Conclusion

This paper proposes a sufficient condition for the global dissipative properties of GRN with external inputs. For the proof of dissipative properties, a large number of documents adopt the method of constructing Lyapunov function. However, this paper adopts the method of matrix measurement, which simplifies the calculation difficulty, and the result is satisfactory. At the same time, error feedback control is given to ensure the stability of the system. And, the limitation of this kind of control is very small, which can bring great convenience in practical application. In addition, the study of the influence of Brown process and Levy process on gene regulatory networks is a challenging problem, which requires more mathematical knowledge about stochastic processes, but it is a very meaningful topic.



---

*Source: 9940229-2021-07-26.xml*