Dynamics of Random Boolean Networks under Fully Asynchronous Stochastic Update Based on Linear Representation

A novel algebraic approach is proposed to study dynamics of asynchronous random Boolean networks where a random number of nodes can be updated at each time step (ARBNs). In this article, the logical equations of ARBNs are converted into the discrete-time linear representation and dynamical behaviors of systems are investigated. We provide a general formula of network transition matrices of ARBNs as well as a necessary and sufficient algebraic criterion to determine whether a group of given states compose an attractor of length in ARBNs. Consequently, algorithms are achieved to find all of the attractors and basins in ARBNs. Examples are showed to demonstrate the feasibility of the proposed scheme.


Introduction
The gene regulatory networks (GRNs) as complex systems present diverse dynamical behaviors. In 1969, random Boolean networks (RBNs) were first introduced by Stuart Kauffman to model gene regulatory networks [1]. Each gene is represented by a node in RBNs with two possible states, i.e. the logical 0 and 1. Generally, the ''1'' value represents the state ''on'' corresponding to a gene that is being transcribed and the ''0'' value represents ''off'' corresponding to a gene that is not being transcribed [2]. A directed edge from one node to another stands for the interaction between genes, the mutual regulation of which is described by an assigned Boolean function. The values of all nodes in a network are updated in parallel and all of them together constitute the state of the whole network. Since dynamics is deterministic and the state space is finite, the states of a particular RBN eventually converge into a series of periodically recurring states, which are called attractors including the fixed points and cycles, starting from an arbitrary initial condition. Attractors and their transient states compose the basins of attraction, which present the dynamical behaviors of systems. RBNs provide a way to formalize large scale complexity in a realistic manner by describing the gene expression with ON/OFF and reflect the statistical features of real living cells by adjusting the parameters of networks [3]. Recent years, along with the research of networks in molecular biology, chemistry, neurobiology, and economy etc [4][5][6][7][8][9][10][11][12][13][14][15][16], RBNs have been extensively investigated and a wealth of results have been achieved [17][18][19][20][21][22][23][24][25][26][27][28].
As abstract models of physical and biological phenomena, a certain degree of simplification is necessary. Usually, synchrony idealization is adopted as update scheme in the classical RBNs, as Ref. [22], ''The simplest class of Boolean network is synchronous, which means that all elements update their activities at the same moment''. All of elements in RBNs evolve according to a global synchronized clock, which is based on a tacit assumption that update scheme would not affect the essential properties of dynamics. However, some doubts have been cast on the validity of the above assumption when more and more evidences showed the dynamic behaviors of systems under asynchronous update presented considerable divergence comparing with the counterparts under synchronous update [29][30]. ARBNs was firstly proposed in 1997 by Harvey et al. [31] where only one node was randomly chosen for update at each time step. Since then, many results on ARBNs have been presented and some modified versions of ARBNs were also proposed [32][33][34][35][36][37][38]. We think, if ''synchronous update'' is considered as ''idealization'' [31], Harvey's model as the extreme opposite of synchronous models should be another form of ''idealization''. Motivated biologically, since the expression of genes is not an instantaneous process, but the transcription of DNA and the transport of enzymes may take from milliseconds up to a few seconds [36]. So when modeling the above dynamics by means of a discrete form, at each time step, the condition to strictly limit only one node for update is too strong to describe the real situation.
Although RBNs have been much analyzed, but as logical systems, they are difficult to be investigated analytically. Recently, a new matrix theory called the semi-tensor product (STP) was proposed by Cheng et al. [39], whereby logical equations can be represented as a linear discrete-time dynamic equation where x(t)~D | n i~1 x i is the semi-tensor product of logical variables. In contrast with usual transition matrix expression, Eq. (1) contains complete information of the logical equations of RBNs, based on which dynamics of systems, such as attractors and basins, can be analyzed. Furthermore, this form can easily be extended to the Boolean control networks. Certain control problems, such as controllability [40,41], observability [42], and realization [43], etc., can be investigated. So far, the previous studies mainly focus on models under synchronous update. However, ARBNs as nondeterministic systems are considerable different from RBNs. For instance, in Ref. [39], Cheng et al. proposed a method by STP technique to find attractors of RBNs, but, which is based on a fact that RBNs are deterministic systems and attractors are circular. Obviously, this prerequisite can't be satisfied for ARBNs. So, we think it's interesting and meaningful to extend the related research into the field of asynchronous random Boolean networks.
Based on the above discussion, firstly, from the aspect of model selection, most of the previous results on ARBNs are based on the Harvey's model. As mentioned above, we think Boolean networks with a more general asynchronous update schedule are worth investigating. Secondly, from the aspect of methodology selection, despite of the overload of matrix operations, STP technique provides a new perspective to make the analytical analysis of logical dynamics. From the viewpoint of theoretical part, it's a valuable technique to be considered. Therefore, in this article, we focus on dynamics of ARBNs where a random number of nodes can be updated at each time step. Mainly, two parts of work are involved in this article. First of all, we complete an algebraic representation of the studied logical models. At the best of our knowledge, the related result is firstly presented. Network transition matrix takes the essential role for the linear representation of logical dynamics based on STP technique. Actually, in the previous works, Cheng et al. [39] have achieved results of RBNs. Here, we engage in a particular class of ARBNs, from the perspective of update schedule, this special update scheme can involve the case of synchronous update, which allow all of nodes are updated at a certain time step. As a result, a general formula of network transition matrices is achieved, based on which one can calculate all of network transition matrices of a particular ARBN including the counterparts under synchronous update and Harvey's model. Secondly, a necessary and sufficient algebraic criterion is discussed, which determine whether a given group of s states could compose a loose attractor of length s in a specific ARBN. Consequently, a novel approach to detect the attractors and basins of ARBNs is designed, which is completely based on set operations and different from the previous studies.
As we know, to search attractors in RBNs is a NP-hard problem. Recent years, many powerful algorithms [3,17,20,35,45] have been presented to improve the search efficiency and reduce the time complexity. Yet, it should be noted that our main interest in this article is to provide a complete quantitative method and a novel perspective to analyze the dynamics of ARBNs from the view point of theoretical part. And results would be taken as preparations for the further research, e.g. Boolean control network under asynchronous update.
This article is organized as follows. In the section of Methods, some concepts and properties of STP are introduced. In Results and discussion, the linear representation of ARBNs is discussed, based on which the dynamic properties of ARBNs are studied and algorithms to find attractors and basins of ARBNs are presented. Some examples are shown to illustrate the main results. Finally, a concluding remark is given.

Methods
In this section, the semi-tensor product of matrices is briefly introduced. Some concepts and properties related to this article are presented.

Definition 1 ([39])
1) Let X be a row vector of dimension np, and Y be a column vector of dimension p. Then we split X into p equal-size blocks as X 1 , . . . ,X p , which are 1|n rows. Define the STP, denoted by , as 2) Let A[M m|n and B[M p|q . If either n is a factor of p, say nt~p and denote it as A[ t B, or p is a factor of n, say n~pt and denote it as A] t B, then we define the STP of A and B, denoted by C~A D |B, as the following: C consists of m|q blocks as C~(C ij ) and each block is where A i is the i-th row of A and B j is the j-th column of B.
It's easy to verify that for two column vectors X [R m and Y [R n , X D |Y [R mn . Note that STP of matrix can be seen as a generalization of the conventional matrix product, so the properties of matrix product, such as distributive rule, associative rule, etc, still hold.
Some notations in this article are defined as follows 1) d r n denotes the r-th column of the n|n identity matrix I n and D n :~d r n 1ƒrƒn j È É , which is the set of all n columns of I n . An operation is defined as H(d r n )~r f g, furthermore, when Next, we define the swap matrix W ½m,n , let X [R m and Y [R n be two column vectors where W ½m,n is a mn|mn matrix labeled columns by (11,12, . . . ,1n, . . . ,m1,m2, . . . ,mn) and rows by (11,21, . . . ,m1, . . . ,1n,2n, . . . ,mn), the elements in position ((I,J),(i,j)) is , I~i and J~j 0, otherwise & : W ½m,m is briefly denoted by W ½m .
In order to get the matrix expression of logic, the Boolean values should be denoted as vector forms T~1*d 1 2 and F~0*d 2 2 . For a logical function L(P 1 , . . . ,P r ) with arguments P 1 , . . . ,P r , it has been proved in Ref. [39] that L(P 1 , . . . ,P r ) can be represented as L(P 1 , . . . ,P r )~M L P 1 P 2 . . . P r , where matrix M L [M 2|2 r is unique, which is called the structure matrix of logical function L.
More details on STP can be found in Ref. [39]. In this article, the matrix products are assumed to be STP and the symbol is omitted.

Model
In this section, the algebraic representation of asynchronous random Boolean networks is discussed. Without special note, ''asynchrony'' in this article implies m t (m t [ 1, . . . ,N f g ) nodes are randomly selected to be updated at time step t. Here, we don't consider the case m t~0 i.e. none of nodes are updated.
ARBNs with N nodes can be described as follows: Where x i (t) represents the state of node i at time t, which receives inputs from p distinct neighbors and x i j' (t) is the j' th input. The nodes in ARBNs take values from the set 0,1 f g corresponding to two levels of gene expressions. f i is a Boolean logic from 0,1 f g p ? 0,1 f g, which is assigned to node i. (2) can be converted into When the in-degree of node i is less than N, a dummy matrix E d~d2 1,2,1,2 ½ should be applied. It's easy to verify that for any Boolean variable u, one can get E d u~I 2 . Assume there was no input from node q to i, x q (t) can be substituted by E d x q (t).
Multiplying all the equations in (3) as According to Ref [39], when  (4) can be represented as where L H t~P This completes the proof.
One can prove that the mapping D | N i~1 : D N 2 ?D 2 n is a bijective mapping. So Eq. (5) can adequately describe the dynamics of Eq.
(3) at time t, say, L Ht involves all of transition information of system from time t to tz1 when set H t is chosen for update, which is called network transition matrix under update elements H t .
Whether RBNs or ARBNs, Network transition matrix is an essential concept for linear representation of logical dynamics based on STP technique. As deterministic system, there is only one network transition matrix for a particular RBN. However, multiple matrices exist in ARBNs depending on different asynchronous update scheme. For models discussed in this article, there are total 2 N {1 network transition matrices. By means of Theorem 1, when getting the structure matrix M i of each Boolean logic f i for any specific ARBN, all of network transition matrices could be calculated.
Example 1. A logical dynamics is described as where x i (t) represents the state of node i at time t, and the symbols ''^'', ''_'' represent conjunction and disjunction, respectively. By truth tables, one can get _*M d~d2 1,1,1,2 ½ and *M c~d2 1,2,2,2 ½ . The algebraic form of Eq. (6) as Case 1: at time t, when only node 1 is selected for update, Case 2: at time t, when all of three nodes are selected for update, Here, we just show the calculations of network transition matrices under the above two cases. Similarly, all of the network transition matrices can be calculated and results are as follows. And all of network transition matrices compose a set T~L 1 ,L 2 ,L 3 ,L 1,2 ,L 1,3 ,L 2,3 ,L 1,2,3 f g . It should be noted that matrix operations are involved in quantity during the above calculations. Fortunately, the related matrices are typically sparse and many efficient algorithms can be applied to optimize the calculation procedure.

The asynchronous network transition table
In order to present the algorithm to find the attractors and basins of ARBNs, we first provide the following definition.
Definition 2. For a random Boolean network with N nodes under asynchronous stochastic update, T~A 1 ,A 2 , . . . , The asynchronous network transition table, briefly called transition table, is defined as Example 2. Recall Example 1. The transition table of system (7) is as follows I~1   2 3 8 1 2 3 8  1 4 1 4 7 8 7 8  1 4 1 4 5 8 5 8  1 4 1 8 3 4 3 8   1 1 4 8 1 1 4 8  1 3 2 4 7 7 8 8  1 3 2 8 3 3 4 8 ð9Þ It should be noted that the asynchronous network transition table is specially designed for ARBNs. As mentioned above, more than one network transition matrix exists in ARBNs, which is different from RBNs. Therefore, there is a question: how to organize these network transition matrices can provide an efficient way for applications. For the above consideration, the asynchronous network transition table is presented, by means of which algorithms to find attractors and basins of ARBNs can be achieved completely based on set operations, which are shown later in Algorithms 1 and 2.

Fixed Points (point attractors)
Proposition 1. For a given ARBN with N nodes, when C(I) i~i f g, the state x~d i 2 N is a fixed point, where I is the asynchronous network transition table.
Proof. Assume T is the set of network transition matrices. When C(I) i~i f g, we can obtain H(Col(L) i )~i f g,VL[T, i.e. the i-th column of L is equal to d i 2 N , which holds for any network transition matrices in T. And, it's easily verified that Ld i

Attractors (Loose attractors)
Fixed points can be seen as point attractors. In the following, we focus on loose attractors in ARBNs. (Sufficiency) Firstly, we should prove that each state in set C can find a way to reach all of the other states. Assume there exists a state d p 2 N [C and all of states that d p 2 N can't reach compose a non-empty set C'. Since d p 2 N can reach each state in set C{C', so any states in set C{C' also can't reach states in C', i.e. C(I) H(C{C') \H(C')~W, which is a contradiction to Eq. (11). So, each state in C can reach the other states.

Theorem 2. For a given ARBN with N nodes and I is the asynchronous network transition table, a state set of
Secondly, we should prove each state in set C can't reach any states out of set C.
We assume there exist d q The proof is complete. As a consequence, we develop a procedure to find all of attractors of ARBNs based on the above result, that is shown later in Algorithm 1.

Basins of attractors
For nondeterministic systems, a partition can't be found to separate the whole state space of ARBNs into several disjointed sets. As mentioned in Ref. [31], basins of attractors in asynchronous random Boolean networks involve two kinds: the definite basins and possible basins. Then, the overlapping parts of basins of different attractors as common phenomena exist the state space of ARBNs, which mean certain transient states can jump from one basin of attractors into another. More studies related to basins in ARBNs can refer to Refs. [32][33][34].
Generally, by a backstepping manner, one can find the basin of a given attractor C. Starting from the state set of attractor C, we can determine parent states of each state in C, which compose a new set C'. Repeat to find the parent states of each state in C' and add them into C' until there aren't any new states involved. The set C'|C is the basin of attractor C. The above procedure is shown in Algorithm 2.
Example 4. Recall Examples 2. The attractors of system (7) include the fixed points 000 and 111. For the fixed point d 1 8 , by checking the transition table I in (9)

Algorithms for detecting the attractors and basins of ARBNs
Given a specific logical dynamics under asynchronous update, all of the network transition matrices could be calculated based on Theorem 1. Consequently, the asynchronous network transition table can also be achieved. By means of the transition table, we formulate the following algorithms to find attractors and basins of ARBNs.  ,s~1 /* F is a feasible state set and s is the length of attractor. */ for (i~1; iv~2 n ; izz) do /*find the fixed points in ARBNs */ if C(I) i~i f g then /* d i 2 n is a fixed point */ Find the basin of state d i 2 n ; /* the procedure to find basin of a particular attractor refers to Algorithm 2*/ Remove the state d i 2 n and its basin from set F ; end if end for while is the number of elements in set. */ Update F ' by removing states whose successors aren't in F '.
while Nonempty (F ') do Randomly pick up one state d i 2 n from F ' and all reachable states from d i 2 n compose a set C. if card(C)~s and H(C)~C(I) H(C) then /* C is an attractor of length s */ Find the basin of attractor C and remove the corresponding states from sets F and F '.
else remove state d i 2 n from F ' end if end while end while End: Algorithm 1 Since, according to the length, the attractors are found successively from small to large, so it isn't necessary to check the condition of Eq. (11) to determine attractors. Algorithm 2: An algorithm for finding basin of a particular attractor C in ARBNs with n nodes.
Begin: Algorithm 2 C'~C /* C is the set of a particular attractor.
end while End: Algorithm 2 In Ref. [39], Cheng's approach to find the attractors of length s in RBNs needs the computation of L s , where L is the network transition matrix. Moreover, Cheng's algorithm is based on a prerequisite that RBNs are deterministic systems and a state on a cycle of length s will be back afters steps. However, this prerequisite can't hold for ARBNs. Therefore, we propose a novel way to detect attractors and basins of ARBNs by means of the transition table.

Examples
In this section, some examples are presented to illustrate the main results of this article. Firstly, an idealized example is shown as follows, which is a model of phosphorylation/dephosphorylation cycles considered by Gonze et. al [44] and re-investigated in [45].
Example 5. A logical dynamics is described as where x i (t) represents the state of node i at time t. System (12) presents a model of protein-protein interactions. As shown in [45], there are three cycles of length four, one cycle of length two, and two fixed points under synchronous update. Here, we study the attractors and basins of system under asynchronous stochastic update.
The algebraic form of Eq. (12) as where x(t)~D | 4 i~1 x i (t). Then, the structure matrices of logical functions are as follows By Theorem 1, all of network transition matrices can be calculated and the asynchronous network transition table is constructed as follows. By means of Algorithm 1, since C(I) 1~1 f g and C(I) 16~1 6 f g in the transition table I, it's easy to get two fixed points fd 1 16 *1111g and fd 16 16 *0000g. Before finding the attractors of system, we first determine the basins of fixed points by Algorithm 2 to reduce the search scope in state space. For the fixed point C~fd 1 16 g, one can get L {1 (C)~d 1 16 ,d 2 16 ,d 3 16 ,d 5 16 ,d 6 16 ,d 9 16 ,d 11 . So far, we can know the basin of the fixed point d 1 16 includes all of the states in state space except for the other fixed point d 16 16 . Therefore, we can conclude that there aren't any attractors in this system. Similarly, the basin of the fixed point d 16 16 involves all of the states except for d 1 16 . A big overlapping part exists between the basins of two fixed points.
One can verify, under synchronous update, System (12) has two fixed points as fd 1 16 g and fd 16 16 g, three cycles of length four as d 2 16 ,d 3 16 ,d 5 (12) is under asynchronous update, all of former cycles turn into the overlapping part of the basins of two fixed points, say, they can jump from one basin of attractor into another. Furthermore, we also can observe some interesting phenomena. Each state in former cycle d 2 16 ,d 3 16 ,d 5 16 ,d 9

16
È É has a direct link pointing to the fixed point d 1 16 , at the same time, all of in-edges of the fixed point d 1 16 come from this former cycle. The same case occurs between the former cycle d 8 16 ,d 15 16 ,d 14 16 ,d 12 16 È É and the fixed point d 16 16 . As to the third former cycle d 4 16 ,d 7 16 ,d 13 16 ,d 10 16 È É , it plays a role of bridge between the rest former cycles of length four, but it hasn't any direct links pointing to either of fixed points. The related state transfers are depicted in Fig. 2. Former cycle d 6 16 ,d 11 16 È É is more special, each state of which has links pointing to all of states in state space except for itself, which is shown in Fig. 3.
Under synchronous update, all of attractors present similar dynamics, e.g. stability, and it's difficult to determine whether there is difference among them. However, some divergences among these attractors can be observed under asynchronous update. Such as the above example in Fig. 2, two former cycles of length four are close to the fixed points, which have greater possibility to flow into the corresponding fixed points. So, we can reasonably consider that the states of the two cycles would be faster into stable states. Yet, the remaining period-4 cycle is in the intermediate position, statistically, whose location determines it would take more time to converge into stability. As to the former period-2 cycle in Fig. 3, to some extent, its characteristic is more unstable under asynchronous update.
Then, the structure matrices of logical functions are as follows.  Based on the above discussions, one can calculate the asynchronous network transition table of system. The detailed computations are skipped and the results are presented as follows.
It can be verified that only one fixed point d 124 256 *(10000100) exists in system (16) under asynchronous update, which is the same as the one under synchronous update. Fig. 5 and Fig. 6 show the basins of the fixed point under asynchronous and synchronous update, respectively.
Comparing Fig. 5 with Fig. 6, under asynchronous update, the basin of the fixed point contains fewer states but has more complex interconnections. There are 11 transient states in Fig. 5 Fig. 6, there are more states directly flowing to the fixed point, furthermore, any states can reach the fixed point though at most 2 steps. There are no attractors in system (16) under asynchronous update, or from another perspective, the remaining states together construct a complicated and big loose attractor.

Conclusions
In this article, we have presented a new approach to study the dynamics of random Boolean networks under asynchronous stochastic update. By semi-tensor product of matrix, the logic of ARBN is converted into linear representation. A general formula of network transition matrices is presented, by which one can obtain all of network transition matrices of a given logical dynamics whether under synchronous update or under asynchronous update. A necessary and sufficient algebraic criterion has been proved to determine whether a given group of states compose attractor of length s on a particular ARBN. Asynchronous network transition table as a transformed set of network transition matrices is presented. Based on the above results, algorithms to detect the  attractors and basins of ARBNs are achieved. Examples are shown to demonstrate the validity of results. As a preparation, results in this article can be applied in the investigation of Boolean control problems under asynchronous update, which would be the work in the future.