Abstract

Here we develop a novel network control framework that integrates the structural and functional information available for intracellular networks to predict control targets. Formulated in a logical dynamic scheme, our approach drives any initial state to the target state with 100% effectiveness and needs to be applied only transiently for the network to reach and stay in the desired state. We illustrate our methods potential to find intervention targets for cancer treatment and cell differentiation by applying it to a leukemia signaling network and to the network controlling the differentiation of helper T cells. We find that the predicted control targets are effective in a broad dynamic framework. Moreover, several of the predicted interventions are supported by experiments.

Author Summary

Here we present a network control approach that integrates the structural and functional information of the network. We show that stabilizing the expression or activity of a few select components can drive the cell towards a desired fate or away from an undesired fate. We demonstrate our methods effectiveness by applying it to a type of blood cell cancer and to the differentiation of a type of immune cell. Overall, our approach provides new insights into how to control the dynamics of intracellular networks. This is a PLOS Computational Biology Methods paper.

Introduction

An important task of modern molecular and systems biology is to achieve an understanding of the dynamics of the network of macromolecular interactions that underlies the functioning of cells. Practical applications such as stem cell reprogramming [1—3] and the search for new therapeutic targets for diseases [4—6] have also motivated a great interest in the general task of cell fate reprogramming, i.e., controlling the internal state of a cell so that it is driven from an initial state to a final target state (see references [7—13]).

Some of these models only include information on which cell components (e.g. molecules or proteins) interact among each other, i.e., the structure of the underlying interaction network. Other models, known as dynamic models, include the structure of the interaction network and also an equation for each component, which describes how the state of this component changes in time due to the influence of other cell components (e.g. how the concentration of a molecule changes in time due to the reactions the molecule participates in).

There are several reasons for this: (i) combinatorial complexity and the size of the matrices involved makes control theory applicable to small networks only, (ii) linear functions are used for the regulatory functions and it is unclear how the switch-like behavior of many biochemical processes [18, 19] will affect these results, and (iii) the notion of controllability in control theory, i.e. control of the full set of states [14—16] or complete controllability, is different from that in the biological sense, which commonly encompasses only the biologically admissible states[8].

For example, Akutusu, Cheng, Tamura et al. [20—22] extend the framework of control theory to systems with Boolean (switch-like) dynamics and provide some formal results in this setting. In the work of Liu et al. [7] the size limitation of linear control theory is overcome by using a maximal matching approach to identify the minimal number of nodes needed to control a variety of real-world large scale networks. Specifically, for some gene regulatory networks, Liu et al. find that control of roughly 80% of the nodes is needed to fully control the dynamics of these networks [7]. In contrast, experimental work in stem cell reprogramming suggests that for biologically admissible states the number of nodes required for control is drastically lower (five or fewer genes [1—3, 8]). Fiedler, Mochizuki et al. [12, 13] use the concept of the feedback vertex set, a subset of nodes in a directed network whose removal leaves the graph without directed cycles (i.e. without feedback loops). They show that, for a broad class of regulatory functions, controlling any feedback vertex set is enough to guide the dynamics of the system to any target trajectory of the uncontrolled network [12, 13]. As one of their examples, the authors use a signal transduction network with 113 elements and show that the minimal feedback vertex set is composed of only 5 elements.

Yet, as manipulating the activity of even a single intracellular component is a long, difficult, and expensive experimental task, it is crucial to reduce as much as possible the number of nodes that need to be controlled. We hypothesize that integrating network structure with qualitative information on the regulatory functions or on the target states of interest could yield control strategies with a small number of control targets. Qualitative information about the regulatory functions is commonly known (e.g. positive/ negative regulation, cooperativity among regulators, etc. ), and relative qualitative information on the desired/undesired states also exists (e.g. upregulation or downregulation of mRNA levels in a disease state with respect to a healthy state). Thus, we choose a logical dynamic framework as our modeling method [23]. This framework is well suited for modeling intracellular networks: discrete dynamic models have been shown to reproduce the qualitative dynamics of a multitude of cellular systems while requiring only the combinatorial activating or inhibiting nature of the interactions, and not the kinetic details [24—30].

., N, each of which denotes the state of a node (also referred to as node state). The state ON (or 1) commonly refers to above a certain threshold level, while the state OFF (or 0) refers to below the same threshold level. The vector formed by the state of all nodes (01, 02,. . ., 0N) denotes the state of the system (or system/network state). To each node 1/,- one assigns a Boolean function which contains the biological information on how node vi’s inputs influence 0,; these functions are used to evolve in time the state of each element. We use the general asynchronous updating scheme [33, 34, 36] (see Methods), a stochastic scheme which takes into consideration the variety of timescales present in intracellular processes and our incomplete knowledge of the rates of these processes.

The attractors of intracellular networks have been found to be identifiable with different cell fates, cell behaviors, and stable patterns of cell activity [24—30, 39, 40]. In general, the task of finding Boolean network attractors is limited by combinatorial complexity; the size of the state space grows exponentially with the number of nodes N. To address this, we recently proposed an alternative approach to find the attractors of a Boolean network which allowed us to identify the attractors of networks for which a full search of the state space is not feasible [41]. This attractor-finding method is based on identifying certain function-dependent network components, referred to as stable motifs, that must stabilize in a fixed state. A stable motif is defined as a set of nodes and their corresponding states which are such that the nodes form a minimal strongly connected component (e.g. a feedback loop) and their states form a partial fixed point of the Boolean model. (A partial fixed point is a subset of nodes and a respective state for each of these nodes such that updating any node in the subset leaves its state unchanged, regardless of the state of the nodes outside the subset.) It is noteworthy that stable motifs are preserved for other updating schemes because of their dynamical property of being partial fixed points. For more details on the attractor-finding method and the identification of the stable motifs see 81 Text and ref. [41]; for a more formal and mathematical discussion see 82 Text section A or Appendix A of ref. [41].

Repeating this procedure iteratively for each separate stable motif until no new stable motifs are found yields the attractors of the logical model. Formally, the result is a set of network states called quasi-attractors, which capture steady states exactly and are a compressed representation of complex attractors [41]. The network control method we propose here builds on the concept of stable motifs and its relation to (quasi)attractors [41] and takes it much further by connecting stable motifs with a way to identify targets whose manipulation (upregulation or downregulation) ensures the convergence of the system to an attractor of interest. The use of quasi-attractors in our method does not compromise its general applicability, but it does require that certain networks with special types of complex attractors are treated with care when our method is applied. None of the networks we discuss in this work nor any intracellular network models we are aware of fall in this category; for more details see 81 Text, 82 Text, and ref [41].)

fB=AORCORNOTE

f0 = (A AND B) OR D

f0 = (NOT B AND NOT A) OR (D AND NOT A) OR (D AND NOT B) OR NOT E

Stable motifs of a logical (Boolean) network. (a) An example of a logical network indicating the regulatory relationships and the logical update function of each node. (b) The four stable motifs of the logical network in (a) and their corresponding node states. These stable motifs are strongly connected components and partial fixed points of the logical network.

This logical network has four stable motifs (Fig 1(b)): (i) {A = 1, B = 1}, (ii) {A = 0}, (iii) {E = 1}, and (iv) {C = 1, D = 1, E = 0}. Network reduction for each of these stable motif yields four reduced networks, each of which has its own stable motifs, all of which are shown in 81 Fig. For example, the reduced logical network obtained from the first stable motif consists of two nodes (D and E) and has two stable motifs: {E = 1} and {E = 0}. The stable motifs of the remaining three reduced logical networks are, respectively: {E = 1} and {D = 1}; {A = 1, B = 1} and {A = O}; {A = 1} and {A = 0}. Repeating the same network reduction procedure with each of the new stable motifs leads to either a new reduced network or one of four attractors (AL, 1' = 1,. . .,4). The stable motifs obtained from the original network and from each reduced network, and the attractors they lead to are shown in Fig 2. This diagram is a compressed representation of the successive steps of the attractor finding process, which include the original network, the stable motifs of the original network, the reduced networks obtained for each stable motif, the stable motifs of these reduced networks, and so on (see 81 Fig). We refer to such a diagram as a stable motif succession diagram, and we note that it is closely analogous to a cell fate decision diagram. We propose to use this stable motif succession diagram to guide the system to an attractor of interest.

Results

Stable motif control implies network control

The stable motifs’ states are partial fixed points of the logical model, and as such, they act as “points of no return” in the dynamics. Normally, the sequence of stable motifs is chosen autonomously by the system based on the initial conditions and timing. We propose to use our knowledge of the sequence of stable motifs to guide the system to an attractor of interest. We refer to this network control method as stable motif control.

We give the proof of this statement in Lemma 4 and Proposition 6 of 82 Text section B. The number of nodes that need to be controlled can be minimized by removing motifs that do not need to be controlled and by finding a subset of nodes in a motif which can fix the whole motifs state. A step by step description of the stable motif control algorithm is given in Methods. For more details on the motif-remov-al step involved in minimizing the number of control nodes, see 81 Text; for a justification of the steps involved in minimizing the number of control nodes, see 82 Text. S3 Text presents a discussion of the complexity of our methods and mitigation techniques for the most time consuming parts of our methods. Table 1. Intervention targets for each control strategy in the T-LGL leukemia network model.

The control set obtained from the first sequence is then {E = O, A = 1}. For the second sequence, a similar reasoning leads to the same control set,

This set was demonstrated to be an effective control target and set an upper limit in the size of the control set in references [12, 13]. Because removing the feedback verteX set from the network must destroy all cycles, including self-loops, there are two possible minimal feedback verteX sets, {A, B, D, E} and {A, C, D, E}. The number of nodes that need to be controlled in our method is half of the size of the feedback verteX set, a substantial improvement. It should be noted that our method does not guarantee that the resulting control sets are small nor that the control sets are the smallest possible, though our case studies suggest that the resulting control sets tend to be relatively small (between one and five nodes out of more than fifty, see Tables 1 and 2, and ref [45]).

Blocking stable motifs may obstruct specific attractors

Based on the motif-sequence point of view provided by the stable motif succession diagram (Fig 2), we hypothesize that blocking the stable motifs that lead to an attractor will either prevent or make it less likely for the system to reach this attractor. We refer to this network control method as stable motif blocking. The algorithm for the method is given in Methods. The interventions obtained from this method are negations of node states of the target attractor, and as such, have the property of eliminating the intended attractor. However, new Table 2. Intervention targets for each control strategy in the helper T cell network. Th1 stable motif blocking interventions (BTm) Th2 stable motif blocking interventions (BThz) Th17 stable motif blocking interventions (BTh17)

In biological situations (like in our test cases) one commonly has certain molecular markers of cell fate which specify the attractor to a large degree but not at the level of every node. Thus the final state obtained after stable motif blocking may still be consistent with the biological specification of the undesired attractor, making the intervention unsuccessful. We also adopt a stricter definition for a successful intervention: if a longterm but not permanent intervention (i.e. a transient intervention) reduces the number of network states or trajectories that lead to the unwanted attractor, then the intervention is considered to be longterm successful. The best-case scenario would be that the manipulated network has only the desired attractors of the original network (i.e., any but the unwanted attractors), in which case the network will stay in these attractors even if the intervention is stopped.

The first intervention (A = 1) causes the system to have .41 and .42 as its only attractors, and thus, the network is driven towards these attractors and away from the unwanted attractor .43. Furthermore, the network stays in those attractors even after the intervention is stopped, as they are also attractors of the original network, so the intervention is longterm successful. Similarly, the second intervention (E = 0) causes the system to have .42 and .44 as its sole attractors, so it is also a longterm successful intervention. The third intervention (D = 0) only leaves attractor A1 intact, and also gives rise to two new at-tractors. To evaluate if this intervention is longterm successful we compare the probabilities that an arbitrary initial condition ends in A3 with and without the intervention. For the intervened case, we set D = 0 for a long time, then stop the intervention and wait for the network to reach an attractor. We find that the intervention makes it more likely for an arbitrary initial condition to reach .43, so this intervention is not longterm successful. Verification of the methods effectiveness in test cases

This is a broad and increasing domain of application: refs. [24—28] are examples of recent logical models that had experimentally validated predictions, while other examples can be found in the review articles [29, 30].

More specifically, we use our network control framework to predict network control interventions on previously developed logical dynamic models for a leukemia signaling network and for the network controlling the differentiation of helper T cells. We confirm the effectiveness of the predicted stable motif control interventions using dynamic simulations, an independent verification of the result we prove in $2 Text. For the case of stable motif blocking interventions, whose effectiveness is not guaranteed, we use dynamic simulations to test the effectiveness of the predicted interventions.

Cytotoxic T cells are a central part of the immune system’s response to infection. These T cells detect antigens in infected cells and, in response, induce the self-destruction of the infected cells. After fighting infection normal cytotoxic T cells undergo activation-induced cell death (apoptosis), but in T-cell large granular lymphocyte (T-LGL) leukemia cytotoxic T cells avoid cell death and survive, which eventually leads to diseases such as autoimmune disorders.

[28]. This network model consists of 60 nodes and 142 regulatory edges, with the nodes representing genes, proteins, receptors, small molecules, external signals (e.g. Stimuli), or biological functions (e.g. Apoptosis). The T-LGL network is shown in Fig 3 and its logical functions are reproduced in S4 Text. Previous work by Zhang et al. [28] and Saadatpour et al. [46] has shown that in the sustained presence of the external signals IL15, PDGF, and Stimuli (antigen presentation) the system has two attractors: one that recapitulates the survival phenotype and node deregulations seen in T-LGL leukemia, and a second one that corresponds to self-programmed cell death (apoptosis) (see S4 Text for more details about attractor specification).

The result is 7 different stable motifs, each of which is shown in Fig 3 with a different node/ edge color (nodes and edges with multiple colors are part of several stable motifs). The stable motif succession diagram for the T-LGL network is shown in Fig 4. For simplicity we do not include the motifs associated with the node P2 in the succession diagram, as these motifs require the other stable motifs to influence the resulting attractor in the succession diagram.

3 and 4) or T-LGL leukemia (the red stable motif in Figs. 3 and 4). This result agrees with previous studies of T-LGL leukemia, in which it was found that blocking SIP signaling induced apoptosis in leukemic T-LGL cells [28, 47] , a result reproduced by the network model when the state of SIP was set to OFF [41, 46].

The obtained intervention targets for each control strategy are shown in Table 1. Note that some intervention targets may be present in both control strategies (e.g. {SIP 2 OFF} is a target both for apoptosis control and T-LGL attractor blocking). For the case of stable motif blocking one may have the same intervention for blocking two different attractors (e.g. {TBET 2 OFF}), which means that this intervention could block either attractor.

The results of the intervention target validation are summarized in S1 Table. For all the stable motif control interventions we obtain 100% effectiveness in reaching the desired state, both for the case in which the intervention is permanent and for the case in which it is not. This means that all stable motif control interventions are longterm successful, in agreement with our formal results in S2 Text. For example, when fixing S1P 2 OFF the apoptosis attractor is reached for all the initial conditions, indicating that the T-LGL attractor is unreachable. For the case of the stable motif blocking interventions we find that each of them but one (GRB2 2 OFF) is successful in blocking its target attractor or one of its target attractors, though not always with 100% effectiveness. For example, for TBET 2 OFF the apoptosis attractor is reached from 10% of the initial conditions, which is a substantial reduction from the baseline of 62% in the case of no intervention, indicating that this interventions is effective as an apoptosis blocking strategy. We also find that most of the stable motif blocking interventions are effective when the intervention is permanent, but only a few of them are effective when the intervention is temporary.

Thus, we evaluate the success of each single intervention from control sets with more than one node (see S1 Table). We find that one of the 12 single node interventions, Ceramide 2 ON, is 100% effective and longterm successful. Of the remaining 11 single node interventions only a few are successful (Ceramide 2 OFF, SPHK1 2 ON, and PDGFR 2 ON) and/or longterm successful (SPHK1 2 ON and PDGFR 2 ON) but none of them are 100% effective. This result illustrates the benefit of combinatorial interventions over single interventions.

Helper T cells are crucial in the regulation of the immune response in mammals. These T cells release specific cytokines that alter how the immune system responds to external agents, for example, by recruiting specific immune system cells to fight infection, promoting antibody production, or inhibiting the activation and proliferation of other cells. Various subtypes of helper T cells are known, such as Th1, Th2, Th17 and Treg, which are distinguished by a differential expression of specific transcription factors and cytokines.

[48]. This network model has several attractors, which correspond to the known canonical helper T cell subtypes, and also to some hybrid cell types (see [48] and S5 Text). The reachability of each attractor depends on the presence of several external environmental signals (either cytokines or antigen), which are represented as input nodes in the network. For our study we use one of the environmental conditions studied by Naldi et al. (TGFB_e 2 ON, IL2_e 2 ON, and APC 2 ON) [48] because it allows us to explore control targets for all T cell subtypes. The helper T cell differentiation network under the selected environmental conditions consists of 55 nodes and 121 edges and is shown in Fig 5. Its corresponding logical functions are reproduced in S5 Text.

Despite the large size of the succession diagram, a closer look at it gives a simple interpretation: the stable motifs associated with each attractor regulate the characteristic transcription factor of each helper T cell subtype. To check this, we look at the minimal subsets of stable motifs that are sufficient for a sequence to lead to a single differentiated helper T cell subtype (see Fig 6 and S5 Text). We use the stable motif succession diagram and our stable motif control and stable motif blocking strategies to find intervention targets for each helper T cell subtype (see Table 2). To validate the proposed intervention targets we use the same procedure as in the T-LGL leukemia network case (see Methods). We also look at the effect of single node interventions for control sets With more than one node. The results of the intervention targets for the stable motif control, stable motif blocking strategies, and single node interventions are summarized in 82 Table. We find that (i) there is a 100% effectiveness in reaching the desired state for all the stable motif control interventions, (ii) most of the stable motif blocking interventions are successful in blocking their target attractor or one of their target attractors, though not always With 100% effectiveness, and (iii) some single interventions are successful, but none of them are 100% effective.

The control targets transcend the logical modeling framework

The network control approach we propose is formulated in a Boolean framework, which brings up the question of Whether the control targets identified are dependent on the logical modeling scheme. To address this, we translate the studied Boolean network models into ordinary differential equation (ODE) models using the method described by Wittmann et al. [49]. In the ODE models the node state variables 6,. can take values in the range [0, 1]; the differential equa-smooth Hill-type function parameterized by Hill coefficients and threshold parameters, and T,is a timescale parameter. The function f 1. is such that it matches the Boolean function Whenever its inputs 61.1, . . . , a 1.10 are either 0 or 1. Thus, the fixed point attractors of the Boolean model are preserved in the ODE model. Table 3. Experimental support for successful control targets in Tables 1 and 2.

We find that the stable motif control interventions are still 100% effective or very close for both permanent and transient interventions (S3 Table and S4 Table). We also find that the effectiveness of the interventions is mostly unchanged by varying the Hill coefficients (S5 Table), varying the the timescale parameters T,- and thresholds (S6 Table), or fixing the intervened node variables close to but not exactly at the intervention-prescribed values (S7 Table). We finally test single interventions and find that they still underperform combinatorial interventions (S3 Table and S4 Table).

We find that several of the single interventions predicted to be successful in inducing apoptosis of leukemic T cells or in inducing specific T cell types were found to be successful experimentally. The control targets for which experimental support was found, the attractors they lead to, and the references are shown in Table 3. Collectively, these results strongly suggest that the control targets identified by our approach transcend the logical framework.

Discussion

Despite recent advances in network controllability approaches, most of them rely solely on the topology [7, 9, 10, 12, 13] or the dynamics [11, 20—22] of the network. Thus, potentially important effects that depend on the interplay between structure (topology) and function (dynamics), such as combinatorial interactions, are not considered. In this work we proposed a network control approach that combines the structural and functional information of a discrete (logical) dynamic network model to identify control targets. The method builds on the concept of stable motif and its relation to finding attractors [41], and takes it much further by connecting stable motifs with a way to identify targets whose manipulation (upregulation or downregulation) ensures the convergence of the system to an attractor of interest. We illustrated our methods potential to find intervention targets for cancer treatment and cell differentiation by applying it to network models of T-LGL leukemia and helper T cell differentiation.

For example, stable motif control interventions are guaranteed to drive an initial state to the target attractor state with 100% effectiveness, regardless of the initial state, a general result which we prove in 82 Text and corroborate in our test cases (see 81 Table and 82 Table). They are also longterm successful, meaning that the intervention only needs to be applied transiently for the network to reach and stay in the desired state, a general result which we also verify in our test cases (see 81 Table and 82 Table). We attribute these properties to the use of the natural (autonomous) dynamics of the network to control its dynamics.

As shown in 81 Table and 82 Table, only one single-node intervention (namely, Ceramide 2 ON in the T-LGL leukemia network) was able to match the 100% effectiveness of the multi-target interventions. This agrees with recent clinical studies on the advantages of combinatorial over single target interventions [50—52]. Finally, the stable motif control interventions for our case studies target only a few nodes (between one and five out of more than fifty), which matches what is expected from stem cell reprogramming experiments [1—3, 8].

Indeed, our control approach is applicable to any dynamic process that can be captured qualitatively by a Boolean dynamic network model such as ecological community dynamics [53], social dynamics [54, 55] , or disease spreading [56, 57]. The validity of the control targets on the translated ODE models of our two case studies and the experimental support found for several of these targets demonstrates the broader, potentially model-independent reach of our method. Further work is needed to address exactly how to extend the concept of stable motif and our network control approach to continuous models; formalizing our framework to admit an arbitrary number of discrete states and other updating schemes may prove a valuable step in this direction. Taken together, our results provide a novel framework for the control of the dynamics of intracellular networks that combines realistically obtainable structural and functional information of the network of interest. As such, we expect this framework to be significant to a variety of practical applications and to also provide a new avenue to better understand how the complex behaviors of cells in living organisms emerges from the underlying network of biochemical interactions.

Methods

Computational methods

The source code of a Java project that allows the user to reproduce the stable motif succession diagrams and control sets for the test cases analyzed is also freely available on GitHub under the examples folder of the StableMotifs Java library. The generation of the ODE model from the logical model was done using the MATLAB implementation of the method of Wittman et al. [49, 58]; the numerical integration of the ODE models was performed using MATLAB’s ode45 function (see S6 Text for more details). The networks in all figures were created using the yEd graph editor (http://WWW.yworks.com/).

General asynchronous updating scheme

In the general asynchronous scheme, the state of the nodes is updated at discrete time steps starting from an initial condition at t: 0. At every time step, one of the variables is chosen randomly (uniformly) and is updated using its respective function and the state of its regulators at the previous time step While the rest of the variables retain their state. In this way, every possible update order is allowed, and thus, all relative timescales of the processes involved are sampled.

Stable motif control algorithm

For an attractor of interest A, the steps of the stable motif network control method are the following: - Step 1: Identify the sequences of stable motifs that lead to A. These can be obtained from the stable motif succession diagram (see Fig 2) by choosing the attractor of interest in the rightmost part and selecting all of the attractor’s predecessors in the succession diagram. - Step 2: Shorten each sequence 8 by identifying the minimum number of motifs in 8 required for reaching A and removing the remaining motifs from the sequence. This minimum number of motifs can be identified from the stable motif succession diagram (Fig 2); they are the motifs after Which all consequent motif choices lead to the same attractor A.

., aml), find the subsets of stable motif s states 0 2 {Mi}, M,- Q M that, When fixed in the logical model, are enough to force the state of every node in the motif into M. At worst, there Will only be one subset, Which Will equal the Whole stable motifs state M. If any of these subsets is fully contained in another subset, remove the larger of the subsets. In each stable motif sequence 8 = (M 1,. . ., M L), substitute every stable motif MJ- With the subsets of the stable motifs states obtained, that is, S = (01,. . ., 0L).

., OL) create a set of states C by choosing one of the subsets of stable motifs states M kj in each Oj and taking their union, that is, C = MkIU' - -UMkL, M kj E Oj. The network control set for attractor A is the set of node states C A 2 {Ci} obtained from all possible combinations of subsets of stable motifs states M kj’s for every sequence 8. To avoid any redundancy, we additionally prune C A of duplicates and remove each set of node states C,- Which is a superset of any of the other sets of node states Cj (i.e. C Ci). For a pseudocode of each step of the stable motif control algorithm see S7 Text.

Stable motif blocking algorithm

Given an attractor A one is interested in obstructing, the steps to identify potential interventions are the following: - Step 1: Identify the sequences of stable motifs that lead to A. This step is the same as the first step in the stable motif control algorithm, and can be obtained from the stable motif succession diagram (Fig 2). - Step 2: Take each stable motifs state MI in the sequences obtained in the previous step. Create a new set M A With all of these stable motif states, M A 2 {Mi}.

Create a new set BA With the negation of each node state, BA 2 {61-}. The node states in BA and any combination of them are identified as potential interventions to block attractor A. For a pseudocode of each step of the stable motif blocking algorithm see S7 Text.

Intervention target validation

We find that, for our test cases, temporal evolution for 10,000 time steps ensures reaching an attractor from any initial condition considered with stable motif control intervention or without an intervention; to be safe, we choose to evolve for 50,000 time steps in all cases. We repeat this for a large number of initial conditions (100,000) and calculate the probability of reaching each attractor from an arbitrary (uniformly chosen) initial condition. We also look at the probability of reaching each attractor when the intervention is not permanent (i.e. it is transient), that is, we fix the prescribed node states for a large number of time steps, then stop fixing these states and wait for another large number of time steps for the system to reach an attractor. For our test cases, we find that using 10,000 time steps for each evolution stage (with and then without prescribed node states) is enough to preserve the first three digits of the estimated probabilities p Am of reaching the attractor of interest, consistent with what is expected from the standard deviation of the estimated probability p Am. To be safe, we choose to evolve for 50,000 time steps for each evolution stage.

For our test cases, we find that 100,000 initial conditions are enough to estimate the probabilities p Am of reaching the attractor of interest with an error (standard deviation of the estimated probability p Am) of 3' 10'3 [p Attr(1—p Attr)] 1/ 2. Equivalently, if p Am is expressed as a percentage (which we denote as %pAttr for clarity), the error in it is estimated as 3- 10‘3[%pA,,.(100%—%pA,,.)] “2% (e.g. 0.03% for a %pA,,, of 1%, and 0.15% for a %pA,,. of 50%). The number of time steps we use is enough to show no changes in p Am beyond what is expected from the standard deviation of the estimated probability p Am, and is also found to be enough for the initial conditions to reach the attractors when no interventions are applied.

Supporting Information

Details and examples of the attractor finding method and stable motif control algorithm. (PDF) 82 Text. Mathematical foundations of the attractor-finding method and of the stable motif control approach. (PDF) S3 Text. Time complexity and mitigation techniques for the attractor-finding method and the stable motif control approach. (PDF) S4 Text. Logical rules and classification of attractors in the T -LGL leukemia network SS Text. Logical rules, classification of attractors, and analysis of the stable motif succession diagram in the helper T cell differentiation network model. (PDF) S6 Text. Translating the logical network models into ordinary differential equation models, and intervention target validation for the ordinary differential equation models. (PDF) S7 Text. Pseudocode for the stable motif control algorithm and the stable motif blocking

List of references that appear in the supporting information files. (PDF)

Read from left to right, the figure shows the logical network in Fig 1, the stable motifs of this logical network, the simplified networks obtained from tracing the downstream effect of each of the original logical network’s stable motifs, and the stable motifs obtained from these simplified networks. Nodes are colored based on their respective node state: gray for 0, black for 1, and white for nodes whose state is not yet determined. Each large arrow has an associated stable motif sharing the arrow’s color. These large arrows stand for the use of a network reduction technique on the network they start from by tracing the downstream effect of their associated stable motifs on this network. (PDF)

The logical function of each example node is shown above its expanded network representation. Nodes are colored White if they denote normal nodes or complementary node (complementary nodes have a bar above their name, While normal nodes do not), and colored black if they denote composite nodes. For more details see 81 Text and 82 Text. (a) Expanded network representation for normal node C, complementary node C , and their inputs. (b) Expanded network representation for normal node B, complementary node B, and their inputs. (PDF)

The figure shows (a) a two node Boolean network whose logical functions are given by an XOR function, (b) the network’s state transition graph, i.e., all combinations of network states and the allowed transitions between them under the general asynchronous updating scheme, and (c) the network’s stable motif succession diagram. This Boolean network is the simplest example (up to a relabeling of node states) of so-called unstable oscillations. Unstable oscillations refer to a subset of nodes whose node states oscillate in an attractor while their node states are fixed in a different attractor, even though both attractors are the same except for the state of this subset of nodes. In the example Boolean network shown in this figure, we have the states of nodes A and B oscillate tors displaying this behavior; for more details see 81 Text and 82 Text. (PDF)

The figure shows (a) a three node Boolean network that displays incomplete oscillations, (b) the sub-state-space of nodes A and B in the network’s state transition graph (i.e., all combinations of network states and the allowed transitions between them) under the general asynchronous updating scheme, and (c) the network’s stable motif succession diagram. Incomplete oscillations refer to a subset of nodes Whose node states oscillate in an attractor but do not Visit all possible states of their sub-state-space in the attractor. In the example Boolean network shown in this figure, we have the states of nodes A and B oscillate between three subnetwork states {(A = 1, B = O), (A = O, tractors displaying this behavior; for more details see 81 Text and 82 Text. (PDF)

Read from left to right, the figure shows the stable motifs of the logical network in Fig 1, the expanded network representation of the stable motifs (from Which stable motifs are formally defined), and the terms of the logical function associated to each stable motif. For more details on the expanded network representation see 81 Text and 82 Text. (PDF)

Validation of the intervention targets in Table 1 and single interventions from control sets with more than one node in Table 1 for the T -LGL leukemia network model. The relative apoptosis % change is defined as (Apoptosis %—Normal apoptosis %)/ (Normal ap-optosis %), Where Normal apoptosis % = 62.1% is the percentage of initial conditions that go to apoptosis When no intervention is applied. Interventions marked With T appear in more than one control strategy or target attractor in Table 1. The percentages are significant in the digits shown and have an estimated absolute error (standard deviation of the mean) of 3'10_3[%

Validation of the intervention targets in Table 2 and single interventions from control sets with more than one node in Table 2 for the helper T cell network. The relative attractor % change is defined as (attractor %—normal attractor %)/ (normal attractor %), Where the normal attractor % is the percentage of initial conditions that go to the attractor of interest When no intervention is applied. The normal attractor percentages are 48.6%, 47.5%, 1.3%, and 2.6% for the Th1, Th2, Th17, and Treg helper T cell subtypes, respectively. Interventions marked With T appear in more than one control strategy or target attractor in Table 2. The percentages are significant in the digits shown and have an estimated absolute error (standard deviation of the mean) of 3-10—3[%pAttr(1OO%—%pAttr)] “2 %, Where %pAttr is the percentage shown (e.g. 0.03% for a %pAttr of 1%, and 0.15% for a %pAttr of 50%).

Validation of the intervention targets in Table 1 for the T -LGL leukemia differential equation network model and single interventions from control sets with more than one node in Table 1 for the T -LGL leukemia differential equation network model. The relative apoptosis % change is defined as (Apoptosis %—Normal apoptosis %)/ (Normal apoptosis %), Where Normal apoptosis % = 54.7% is the percentage of initial conditions that go to apoptosis When no intervention is applied. Interventions marked With T appear in more than one control strategy or target attractor in Table 1. The percentages are significant in the digits shown and have an estimated absolute error (standard deviation of the mean) of 6- 10_3[%pAttr(100%—% pAttr)]1/2 %, Where %pAm is the percentage shown (e.g. 0.06% for a %pAm of 1%, and 0.3% for a %PAttr 0f 50%).

Validation of the stable motif control intervention targets in Table 2 for the helper T cell differential equation network model. The relative attractor % change is defined as (attractor %—normal attractor %)/ (normal attractor %), Where the normal attractor % is the percentage of initial conditions that go to the attractor of interest When no intervention is applied. The normal attractor percentages are 50.0%, 45.4%, 2.8%, and 1.8% for the Th1, Th2, Th17, and Treg helper T cell subtypes, respectively. Interventions marked With T appear in more than one control strategy or target attractor in Table 2. The percentages are significant in the digits shown and have an estimated absolute error (standard deviation of the mean) of 6'10—3[%PAttr(100%_%PAttr)]1/2 %, Where %pAm is the percentage shown (e.g. 0.06% for a % pAttr of 1%, and 0.3% for a %pAttr of 50%).

Validation of some stable motif control intervention targets in Table 1 for different Hill coefficients (n) in the T -LGL leukemia differential equation network model. The percentages are significant in the digits shown and have an estimated absolute error (standard deviation of the mean) of 6-10_3[%pAttr(100%—%pAttr)]1/2 %, Where %pAm is the percentage shown (e.g. 0.06% for a %pAttr of 1%, and 0.3% for a %pAttr of 50%).

Validation of some stable motif control intervention targets in Table 1 when fixing the intervened node variables close to the intervention-prescribed values in the T -LGL leukemia differential equation network model. If the intervention is 0 (1), the node variable is fixed at 0.1 (0.9), 0.8 (0.2), 0.7 (0.3), or 0.6 (0.4). The percentages are significant in the digits shown and have an estimated absolute error (standard deviation of the mean) of 6- 10_3[%

Acknowledgments

We would like to thank Steven N. Steinway for fruitful discussions.

Author Contributions

Performed the experiments: IGTZ. Analyzed the data: IGTZ RA. Contributed reagents/materials/ analysis tools: IGTZ. Wrote the paper: IGTZ RA.