Volume 14 Supplement 10
Selected articles from the 10th International Workshop on Computational Systems Biology (WCSB) 2013: Bioinformatics
Computing preimages of Boolean networks
 Johannes Georg Klotz^{1}Email author,
 Martin Bossert^{1} and
 Steffen Schober^{1}
DOI: 10.1186/1471210514S10S4
© Klotz et al; licensee BioMed Central Ltd. 2013
Published: 12 August 2013
Abstract
In this paper we present an algorithm based on the sumproduct algorithm that finds elements in the preimage of a feedforward Boolean networks given an output of the network. Our probabilistic method runs in linear time with respect to the number of nodes in the network. We evaluate our algorithm for randomly constructed Boolean networks and a regulatory network of Escherichia coli and found that it gives a valid solution in most cases.
Introduction
In systems and computational biology Boolean networks (BN) are widely used to model regulative dependencies of organisms [1, 2]. We consider networks, which map a set of environmental conditions to the presence of proteins and finally to actual chemical reactions, which are often modeled as fluxes of a fluxbalance analysis [3]. Hence, these networks are used to make in silico predictions of behavior of organisms in a certain environment [4].
In this paper we address the inverse problem, i.e., we want to predict environmental conditions that allow certain reactions to take place, and others not. Hence, in general, we need to find a set of possible inputs that lead to a given output. This so called predecessor problem or preimage problem has been addressed by Wuensche in [5] and has been shown to NPhard in general [6], which makes it infeasible to solve it for large networks. In [7] an algorithm with reduced complexity for BNs with canalizing Boolean functions has been introduced. However, the problem is infeasible under certain conditions. Both algorithms are designed to find the whole set of preimages, i.e., all inputs to the BN with lead to a certain, desired, output.
In some applications, knowledge of the whole preimage set is not important, merely it can be sufficient to know a subset of the preimage set. Here, we propose a probabilistic algorithm, which solves this problem in linear time with respect to the number of nodes in the network, based on a variation of the well known SumProduct Algorithm (SPA) [8], which is used for a variety of tasks, including decoding error correction codes in communication engineering [9].
Methods
Boolean networks and main idea
The network itself consists of n nodes, and a set of directed edges connecting these nodes. Each node i has a certain state, which can be either zero or one, represented by a variable x_{ i }. Its value is determined by evaluating a Boolean function (BF) f_{ i }. Further, lets define the set ñ(f_{ j }) as the incoming nodes of node j. For example in Figure 1, ñ(f_{5}) = {1, 3}. The BF f_{ j } is a function mapping k_{ j } = ñ(f_{ j }) values of {0, 1}^{ k } to {0, 1}, where k is also called the indegree of node j. The number of edges emerging from a node is called outdegree.
as the wellknown SPA can be used to compute the marginals efficiently. If the approximation is good enough sampling out the product of the marginals will yield an element in Ω_{ y } with reasonable probability.
Proposed algorithm
In this section we will first discuss the basic principles of factor graphs and the SPA. Then we will describe the BN as factor graph and will formulate the actual algorithm to find the marginals. Finally, the sampling is described.
Factor graphs and sumproduct algorithm
where X_{ j } is the subset of [n] containing the argument of h_{ j } . We can then define a factor graph [8] as a bipartite graph consisting of n nodes representing variables {x_{1}, . . . , x_{ n }} (variable nodes) and of m nodes representing functions {h_{1}, . . . h_{ m }} (function node). Edges exist between a function node and a variable node if and only if x_{ i } is an input to function h_{ j } .
In general the computation of the g_{ i } is difficult, but due to the factorization of g the task can be efficiently solved using the the so called SumProduct Algorithm (SPA) [8]. The algorithm iteratively passes messages between the nodes of the graph. At each iteration the messages µ are sent from the function nodes to the variable nodes, containing the corresponding marginal function of the local function. These messages are computed as follows [8]:
Function to variable node
where n(i) give the set of neighboring nodes of node i.
At the variable nodes, these messages are then combined to a marginal function λ and sent back to the function nodes [8]:
Variable to function node
The Boolean network as factor graph
We apply the concept of factor graphs to BNs. Each node in the network represents one variable x_{ i } ∈ {0, 1}, i ∈ [n] of the factor graph, hence we have n variable nodes. Each BF f_{ j } of the BN $\left(j\in \left[n\right]\backslash \mathbb{I}\right)$ is a function node and is connected to the node j and the incoming nodes ñ(f_{ j }). Lets to define ${\tilde{X}}_{j}$ as the variables of the incoming nodes of node j, i.e. the argument of the BN f_{ j }. Further, we define ${\tilde{X}}_{j}^{\left(i\right)}$ as ${\tilde{X}}_{j}$ without the node i.
This problem is an instance of the problem described in Section Factor Graphs and SumProduct Algorithm, hence we apply the same methods here.
Update rule: function to variable node
Lets define ñ(f_{ j }) as the set of indices of the input nodes of the BF f_{ j }.
where λ_{ l } is the probability distribution of variable node l and is defined in Eq. 3.
Update rule: variable to function node
where ${S}_{j}$ is the set of all function nodes, which have node j as input.
Finding the input distributions
A scheme of the algorithm is given in Algorithm 1.
The probability distribution of each node j ∈ [n] at iteration t is given as ${L}_{j}^{\left(i\right)}$ and are initialized with ${L}_{j}^{\left(0\right)}=0$, which is equivalent to the uniform distribution. Then we set the LLRs for the outnodes to either −∞ or +∞ depending on the desired output y of the BN. At each iteration the algorithm can be split in two steps. The first step iterates over all function nodes $\left(j\in \left[n\right]\backslash \mathbb{I}\right)$ and all input variables i ∈ ñ(f_{ j } ) calculating the LLR ${L}_{j\to i}^{\left(t\right)}$ using Eq. (2) and Eq. (4).
In the second step we update all variablesnodes, where the LLRs L_{ j } represents the distributions λ_{ j } and, hence the product of Eq. 3 becomes a summation. Please note, that the LLR of the previous iteration is also added to the sum, in order to prevent rapid changes of the distributions.
After performing a certain number of iterations t_{ max }, the desired marginal distributions of the input variables are found.
Algorithm 1
Initialize ${L}_{j}^{\left(0\right)}=0$ for all nodes
Set the desired LLRs of the outnodes, i.e., ${L}_{j}^{\left(0\right)}$ is either −∞ or +∞, for all outnodes $j\in O$.
t = 0
repeat
t=t+1
for each noninnode $\left(j\in \left[n\right]\backslash I\right)$do
for each input variable i ∈ ñ(f_{ j }) do
calculate L_{ j→i } using Eq. (2) and Eq. (4)
end for
end for
for each nonoutnode v do
${L}_{j}^{\left(t\right)}={L}_{j}^{\left(t1\right)}+{\sum}_{l\in {S}_{j}}{L}_{l\to j}^{\left(t\right)}$
end for
until maximum number of iterations reached
Sampling
The sampling part of our approach is straightforward. Using the marginal distributions ${L}_{j}^{\left({t}_{max}\right)}$, $j\in \mathbb{I}$ we randomly draw vectors x and check if they fulfill y = f(x). If so, they are added to the set ${\tilde{\text{\Omega}}}_{y}$. This procedure is repeated for a certain number of samples.
Simulation results and discussion
We tested our algorithm with randomly generated networks and the regulatory network of Escherichia coli (Ecoli) [2]. The random networks consist of 2400 nodes with N = 200 and M = 1200. We have chosen the BFs from:
· all functions with k ≤ 15 (Type A)
unate, i.e. locally monotone, functions with k ≤ 15 (Type B)
One can see, that for t_{ max } ≥ 14 there is almost no improvement in the similarity. This number is equal to two times the number of nodes between input and output, i.e., it seems to be sufficient that the messages travel once through the network and back. Thus, the following simulations have been perform setting t_{ max } = 14.
Simulation results for different networks
network  num of samples  solved  valid  Unique 

Type A  1000  89%  608.81  4.43 
Type B  1000  95.9%  270.74  68.60 
Ecoli  1000  98.6%  193.3  193.3 
One can see from the results, that in general for most networks and y s at least one preimage can be found. It is worth mentioning, that for the Ecoli network every sampled solution was unique. This is due to the fact, that there exist a few inputs, who completely determine the output. The other input variables have then no influence and hence a marginal distribution of 0.5. Further, the results for the network of type B are much better than for type A. It seems that the marginal distributions for unate functions give better estimation of the actual distribution than the marginal distributions for nonunate functions.
Conclusions
In this work, we proposed a probabilistic algorithm to address the preimage problem of Boolean networks. This is of interest when designing experiments, in which certain regulators are supposed to be in a specific state. Performing a series of simulations with Random networks we showed, that the algorithm works not only for unate functions, of which most biologically motivated networks consist, but for any kind of Boolean functions. By replacing the fixed output values of the network by probabilities one can simply apply the algorithm to networks, whose designated output is described by probability distributions. Further, the algorithm may be easily adjusted to work on stochastic, e.g. Bayesian, networks, where the nodes contain only transition probabilities instead of Boolean function. Therefore, it is needed to adapt the update rules accordingly. It remains an open question, which influence topographic properties, such as number of layers and number of nodes in these layers, have to the performance of the proposed algorithms, since we only investigated networks which are similar to the regulatory network of Ecoli.
List of abbreviations
 Ecoli :

Escherichia coli
 BF:

Boolean Function
 BN:

Boolean Network
 Eq:

Equation
 LLR:

LogLikelihood Ratio
 SPA:

SumProduct Algorithm
Declarations
Acknowledgements
The authors would like to thank Shrief Rizkalla for implementing parts of the simulation.
Declarations
This work was supported by the German research council " Deutsche Forschungsgemeinschaft" (DFG) under Grant Bo 867/252 and the Ulm University.
This article has been published as part of BMC Bioinformatics Volume 14 Supplement 10, 2013: Selected articles from the 10th International Workshop on Computational Systems Biology (WCSB) 2013: Bioinformatics. The full contents of the supplement are available online at http://www.biomedcentral.com/bmcbioinformatics/supplements/14/S10
Authors’ Affiliations
References
 Kauffman SA, Peterson C, Samuelsson B, Troein C: Random Boolean network models and the yeast transcriptional network. Proceedings of the National Academy of Sciences of the United States of America. 2003, 100 (25): 1479614799. 10.1073/pnas.2036429100.PubMed CentralView ArticlePubMedGoogle Scholar
 Covert MW, Knight EM, Reed JL, Herrgard MJ, Palsson BO: Integrating highthroughput and computational data elucidates bacterial networks. Nature. 2004, 429 (6987): 9296. 10.1038/nature02456.View ArticlePubMedGoogle Scholar
 Feist AM, Henry CS, Reed JL, Krummenacker M, Joyce AR, Karp PD, Broadbelt LJ, Hatzimanikatis V, Palsson BO: A genomescale metabolic reconstruction for Escherichia coli K12 MG1655 that accounts for 1260 ORFs and thermodynamic information. Molecular Systems Biology. 2007, 3: 121PubMed CentralView ArticlePubMedGoogle Scholar
 Feuer R, Gottlieb K, Viertel G, Klotz J, Schober S, Bossert M, Sawodny O, Sprenger G, Ederer M: Modelbased analysis of an adaptive evolution experiment with Escherichia coli in a pyruvate limited continuous culture with glycerol. EURASIP Journal on Bioinformatics and Systems Biology. 2012, 2012: 1410.1186/16874153201214.PubMed CentralView ArticlePubMedGoogle Scholar
 Wuensche A: The Ghost in the Machine: Basins of Attraction of Random Boolean Networks. Artificial Life III Proceedings, Santa Fe Institute Studies in the Sciences of Complexity. 1994, AddisonWesleyGoogle Scholar
 Akutsu T, Hayashida M, Zhang SQ, Ching WK, Ng MK: Analyses and Algorithms for Predecessor and Control Problems for Boolean Networks of Bounded Indegree. Information and Media Technologies. 2009, 4 (2): 338349.Google Scholar
 Klotz JG, Schober S, Bossert M: On the Predecessor Problem in Boolean Network Models of Regulatory Networks. International Journal of Computers and Their Applications. 2012, 19 (2): 93100.Google Scholar
 Kschischang F, Frey BJ, Loeliger HA: Factor Graphs and the SumProduct Algorithm. IEEE Transactions on Information Theory. 2001, 47: 498519. 10.1109/18.910572.View ArticleGoogle Scholar
 Gallager RG: LowDensity ParityCheck Codes. 1963, Cambridge: M.I.T. PressGoogle Scholar
 Hagenauer J, Offer E, Papke L: Iterative decoding of binary block and convolutional codes. Information Theory, IEEE Transactions on. 1996, 42 (2): 429445. 10.1109/18.485714.View ArticleGoogle Scholar
Copyright
This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.