BACKGROUND
reverse engineering of biological networks is a key for the understanding of biological systems.
the exact knowledge of interdependencies between proteins in the living cell is crucial for the identification of drug targets for various diseases.
however, due to the complexity of the system a complete picture with detailed knowledge of the behavior about the individual proteins is still in the far future.
nonetheless, the advent of gene perturbation techniques, like rna interference   <cit> , has enabled the possibility to study cellular systems systematically under varying conditions, hence opening new perspectives for network reconstruction methods.
a number of approaches have been proposed in the literature for estimating networks from perturbation effects.
many of these approaches aim at reconstructing a network from directly observable effects.
for example, rung et al.
<cit>  builds a directed disruption graph by drawing an edge , if gene i results in a significant expression change at gene j. wagner  <cit>  uses such disruption networks as a starting point for a further graph-theoretic method, which removes indirect effects  <cit> , hence making the network more parsimonious.
tresch at el.
<cit>  extend this approach by additionally making use of p-values and fold-change directions to make the network more consistent with the observed biological effects.
also bayesian networks have been used to model the statistical dependency between perturbation experiments  <cit> .
for this purpose pearl  <cit>  proposes an idealized model of interventions.
he assumes that once a network node is manipulated, the influence of all parent nodes is eliminated and the local probability distribution of the node becomes a point mass at the target state.
besides for bayesian networks, ideal interventions have also been applied for factor graphs  <cit>  and dependency networks  <cit> .
epistasis analysis offers a possibility for learning from indirect downstream effects.
for example, driessche et al.
<cit>  use expression profiles from single and double knockdowns to partly reconstruct a developmental pathway in d. discoideum via a simple cluster analysis.
also fully quantitative models using differential equation systems have been suggested.
for example, nelander et al.
<cit>  propose a model for predicting combinatorial drug treatment effects in cancer cells.
recently, nested effects models   <cit>  have been proposed as a method, which is specifically designed to learn the signaling flow between perturbed genes from indirect, high-dimensional effects, typically monitored via dna microarrays.
nems use a probabilistic framework to compare a given network hypothesis with the observed nested structure of downstream effects.
perturbing one gene may have an influence on a number of downstream genes, while perturbing others affects a subset of those.
moreover, several of these subsets could be disjoint, i.e.
the knockdown of gene i shows effects, which mostly do not overlap with the effects seen at the knockdown of gene j. nems have been applied successfully to data on immune response in drosophila melanogaster  <cit> , to the transcription factor network in saccharomices cerevisiae  <cit> , to the er-α pathway in human breast cancer cells  <cit> , and to the synthetic lethality interactions network in saccharomicies cerevisiae  <cit> .
the work presented in this paper is designed for a different scenario: we would like to reverse engineer a protein signaling network based on experimentally measured effects on protein expression and activation level after multiple interventions.
these interventions may also be combinatorial  <cit> , i.e.
there is more than one knock-down at a time.
importantly, the set of all perturbations should cover a fraction as large as possible of the network proteins.
effects of all interventions on the network proteins are quantified directly on protein expression and activation level, for instance via reverse phase protein arrays   <cit> .
here, we propose a probabilistic approach called deterministic effects propagation networks , which can estimate the most likely signaling network based on these data.
depns are a special case of bayesian networks, which employ a mixture of purely deterministic and gaussian variables.
while depns and nems have a similar effects propagation scheme, depns differ from nems with respect to the following aspects:
• in nems each node corresponds to one perturbation experiment.
in depns each node corresponds to one single protein, potentially influenced by one or several perturbations.
• nems only work, if the number of perturbations is much smaller than the number of monitored downstream effects.
in depns the opposite is true: depns assume that at least each network node has been perturbed once.
more perturbations are beneficial.
• in nems indirect, high dimensional downstream intervention effects are monitored.
in depns intervention effects are monitored on all other network proteins .
• depns offer the possibility to include combinatorial interventions in addition to single interventions in a natural way.
• depns can use latent network nodes .
this has not been introduced for nems so far.
• depns offer a mechanism for missing value imputation.
this has not been introduced for nems so far as well.
in contrast to other bayesian networks, depns allow for directed graph structures with loops, because they employ a deterministic effects propagation scheme.
our approach was first validated using an extensive simulation study.
afterwards it was employed to reverse engineer the erbb receptor regulated g1/s transition network in hcc <dig> human breast cancer cells, which show de novo resistance against trastuzumab .
in this dataset  <cit>  knock-downs of  <dig> proteins utilizing rnai were performed and effects were measured via rppas.
RESULTS
deterministic effects propagation networks
our basic assumption is that we measure an unknown signaling network of proteins by performing multiple interventions on these proteins.
these perturbations may affect single network proteins at a time or may also be combinatorial, i.e.
they affect several network proteins at one time.
for each network protein in the network we monitor the effect of all interventions carried out on proteins in the network.
note, that each specific intervention not only influences direct targets, but may also cause effects on downstream network proteins.
we explicitly assume that in one experiment all network proteins are unperturbed.
we employ special instances of bayesian networks, which we call deterministic effects propagation networks , to infer the most likely protein network given measurements from multiple interventions under different experimental conditions  or time points.
our implementation of depns can deal with latent network nodes  and with missing data.
the details of our method are described in the section material and methods.
briefly, the idea is that we have an unknown network graph, where each node can have two states: perturbed and unperturbed.
in principle, this could also be extended to three states , but in this paper we only deal with the simpler binary case.
furthermore, attached to each node we have experimental measurements, which for each time point are assumed to come from two gaussian distributions, one for the perturbed and one for the unperturbed case.
the parameters of these distributions are determined in the following way:
<dig>  for each perturbation experiment we compute the expected downstream effects.
this is done by assuming node i to be perturbed whenever i itself or any of its parents are perturbed, i.e.
perturbation effects are propagated in a deterministic fashion from parents to children.
this implies that the network structures learned by depns are always transitively closed graphs: if there is a path from a to b then there is also an edge a → b. it should to mentioned that the perturbations state of a protein can never depend on itself through a loop in our model, i.e.
when the state of a variable is obtained it is not updated further.
however, this behavior could be changed in principle.
<dig>  since we have several perturbation experiments and also one control experiment, we can now determine for each protein measurement, whether the protein was perturbed or not.
<dig>  with this information the data for each protein gets divided into two distributions: one for the perturbed and one for the unperturbed case.
the parameters of these distributions, which are supposed to be gaussian, are estimated either in a maximum likelihood or in a bayesian fashion.
the details are described in the material and methods section.
given the parameters we can then calculate the probability of observing the data for all proteins under a given network hypothesis.
an edge a → b in this network denotes an influence of a on b. this may be interpreted as the direction of the signal flow.
our model does not distinguish between signaling via a transcriptional regulation or via a protein phosphorylation.
the type of the edge is dependent on the type of antibody used .
if the number of nodes in the network is considerably small , we can enumerate all possible network hypothesis and take the one with the highest likelihood.
however, this exhaustive search will be impossible as soon as the network gets larger.
for instance, for  <dig> nodes there are already  <dig> network hypothesis to test.
therefore, here we restrict ourselves to a greedy hill climbing approach, in which we begin with an initial network  and then successively add the edge  increasing the likelihood most.
all possible alternatives are checked.
we do not use any arbitrary starting point for the hill climber, because depns operate in the space of transitively closed graphs only.
therefore, no backtracking/edge removal step can be implemented in an obvious way.
hence, initial networks should be chosen in a sensible way.
simulation studies
graph sampling
we tested the performance of our approach in a simulation study on artificial data, which were created as follows: subgraphs with n =  <dig> and n =  <dig> nodes were cut out randomly from randomly selected kegg  <cit>  signal transduction pathways.
only gene-gene interactions were considered as edges in the graph.
please note that our sampled subgraphs are not random graphs, but randomly selected existing subgraphs of kegg pathways.
for n =  <dig>  <dig> % of the sampled subgraphs were acyclic.
for n =  <dig> these were  <dig> %.
data simulation
single knock-downs of all individual network proteins as well as one experiment without interventions were simulated as follows: for each simulated experiment m replicates were drawn from normal distributions for each network node.
the means of these normal distributions were sampled from n, if the node was expected to be perturbed according to the deterministic effects propagation scheme described above, and from n otherwise.
the variances of the normal distributions were drawn from inv-χ <dig>  the choice of these parameters corresponds to those in section "methods" and is in agreement with our experimental data.
network reconstruction
figures  <dig>   <dig>   <dig> and  <dig> visualize the performance of a network reconstruction with our approach for a varying number of replicates in terms of sensitivity  and specificity .
for each number of replicates  <dig> networks with simulated perturbation data were generated as described before.
as seen from the plots our method achieves a constantly high specificity of >90% with a sensitivity between 70% and 80%, i.e.
the number of replicates has only a minor influence on the quality of the network reconstruction.
the reason for this robust behavior is that due to the effects propagation scheme, we only operate in the space of transitively closed graphs.
the restricted model class can be viewed as a kind of regularization  <cit>  helping to identify a well performing network structure.
additionally, effects propagation increases the chance of getting enough data for both, the perturbed and the unperturbed case, for reliable parameter estimation.
we compared our depn approach to a standard implementation of bayesian networks, which directly models the dependency of measurements : for this purpose we discretized our simulated data .
this was done by setting  <dig>  if the data entry was closer to  <dig>   than to  <dig>  , and  <dig> otherwise.
the conditional distribution at each network node was therefore a binomial with parameter equal to the perturbation probability of the node.
we employed the r package bnlearn to learn the bn network structure via a hill climbing procedure initialized with an empty network.
the bayesian dirichlet score was used to evaluate network structures.
please note, that bnstandard in contrast to depns can only learn network structures without cycles.
however, the vast majority of all original kegg graphs were acyclic here .
since our simulated data contains only single knock-outs, we also compared our depn approach to nested effects models : for this purpose we discretized our simulated data as before.
the likelihood model introduced in  <cit>  was employed with type i and type ii error rates  <dig> , which are expected by the known variances and means of the perturbed and unperturbed distributions.
the same greedy hillclimber as for depns was employed for network structure search.
figures  <dig>   <dig>   <dig> and  <dig> show a similar performance of bnstandard compared to depns in terms of specificity, but a significantly worse sensitivity.
the performance on acyclic graphs was essentially the same as on all graphs .
in terms of computation times bnstandard is clearly faster than the depn approach .
future optimizations of our depns implementation may reduce this gap.
nevertheless, depns are  only suitable for small network reconstructions.
this is, however, in the light of the large experimental effort for multiple interventions and determination of their respective effects on network proteins not really a practical limitation.
in summary the comparison clearly demonstrates that a standard bayesian network approach, which models only conditional dependencies between measurements while ignoring the perturbation signal flow in the network, is not suitable for our kind of data.
compared to nems, depns show a significantly higher specificity and a slightly reduced sensitivity .
in terms of computation times, nems are clearly faster than depns .
this is expected, because nems can be formulated in a much more efficient way  <cit> .
since nems only employ single gene perturbations, they allow to perform deterministic effects propagation in an implicit rather than explicit way and do not have to estimate distribution parameters like depns.
apart from that nems and depns have a lot of similarities, but also major differences, which we highlight in section "related work".
altogether, both comparisons demonstrate that depns seem to be well suited for the low dimensional knock-down effects data, which we are interested in here.
missing values
we next investigated the effect of an increasing number of missing values in a dataset with  <dig> replicate measurements in networks with n =  <dig> and n =  <dig> nodes, i.e.
<dig> ×  <dig> ×  <dig> =  <dig> and  <dig> ×  <dig> ×  <dig> =  <dig> data points in total .
the positions of missing values were chosen uniform randomly in the complete data matrix.
as seen from figures  <dig>   <dig>   <dig>   <dig> our method behaves very robust against replacing missing values with estimated posterior distribution mode values, even with  <dig> out of  <dig>  and  <dig> out of  <dig>  missing values in total.
again, this robust behavior can be attributed to the same reasons as described before: a restricted model class together with the deterministic effects propagation scheme.
latent nodes
we also monitored the effect of having latent nodes in the network.
again, data with  <dig> replicates was simulated in networks with n =  <dig> and n =  <dig> nodes.
the results shown in figures  <dig>   <dig>   <dig>   <dig> reveal a constantly high specificity of >90%, whereas the sensitivity drops with  <dig> latent node from ~80% to ~60% with n =  <dig> and ~80% to ~70% with n =  <dig>  respectively.
with  <dig> latent nodes we still achieve a sensitivity of ~40% with n =  <dig> and > 60% with n =  <dig>  with  <dig> latent out of  <dig> nodes we have ~20% sensitivity, and with  <dig> latent out of  <dig> nodes we have ~40% sensitivity.
we have to remind here that only outgoing edges from latent nodes can be estimated , hence the result shows an expected behavior.
note that on average with each latent node the number of unobservable edges should increase by roughly 1/ <dig> = ~17% with n =  <dig> and by 10% with n =  <dig> .
we believe that a positive point of our method is that still the specificity is rather high, meaning we do not get too many false positive edges.
reconstruction of the erbb receptor-regulated g1/s transition network
rnai and reverse phase protein arrays
sixteen proteins playing a role in the g1/s transition of human cells were taken.
among these proteins there were the receptors erbb <dig>  erbb <dig> and erbb <dig>  several signaling intermediates  and cell cycle proteins .
finally, we included the retinoblastoma protein prb <dig>  which regulates the g1/s transition.
interactions between these proteins known from the literature are shown in figure  <dig> and additional file  <dig>
sixteen knock-downs  with chemically synthesized sirnas for proteins in our network and one experiment with mock transfected cells were conducted.
after each sirna transfection cells were stimulated with egf for  <dig> hours.
we measured the protein expression of signaling intermediates  for each knockdown sample using reverse phase protein arrays  <cit>  before and after egf stimulation with  <dig> technical and  <dig> biological replicates.
these data were integrated using quantile normalization  <cit> .
further details of the experimental setup and a list of performed interventions can be found in the material and methods part of this paper.
additional information can be found in our earlier publication  <cit> .
a heatmap of the data matrix is available in the supplement of this paper .
network reconstruction with depns
we used depns to learn a network between the sixteen proteins purely from data.
since for  <dig> out of the  <dig> proteins in our network no measurements were available due to the lack of antibodies, they were marked as latent nodes .
furthermore, for perk1/ <dig> we had only  <dig> instead of  <dig> biological replicates for cyclin e <dig> and cyclin d <dig> knock-downs.
since our method in the current implementation requires the same amount of data for each experiment, measurements for the third biological replicate were treated as missing values and imputed during network learning using the em algorithm approach described in the materials and methods part of this paper.
network structure learning was performed via greedy hill climbing.
in order to retrieve a measure of confidence for the individual edges learned by our method, we performed a nonparametric stratified bootstrap  <cit> .
that means we sampled  <dig> times with replacement from each replicate group within each condition .
on each bootstrap sample  we estimated a network using depns.
after the bootstrap we counted the relative frequency of each edge in all  <dig> depn reconstructions.
the relative frequency is an approximate probability for each edge .
only edges with a probability >  <dig>  were considered as being reliable.
the result shown in figure  <dig> reveals that our methods correctly recovered 20/ <dig> interactions from the literature network with probability >  <dig> , although we had  <dig> latent nodes.
here we also counted edges as correct, if they correspond to indirect interactions in the original network.
this was valid, since depns learn transitively closed graphs.
the edge erbb <dig> → erbb <dig> could be explained by the fact that erbb <dig> and erbb <dig> together form a heterodimer activating the mapk and pi3k pathway  <cit> .
likewise, the edge cdk <dig> → cdk <dig> was there, because both proteins form a complex  <cit> .
the complex formation of cdk <dig> and cdk <dig> could also explain the edge cdk <dig> → p <dig>  <cit> .
it was found that akt <dig> increases the expression of human igf1r  <cit> , which explained the edge pakt <dig> → perk1/ <dig> as an indirect activation of perk1/ <dig> by akt <dig> through igf1r.
in conclusion there were only three edges, namely cdk <dig> → erbb <dig>  pakt <dig> → erbb <dig> and myc → erbb <dig> for which we could not find any explanation in the literature.
we additionally performed a network inference starting from the literature network as an initial network for the greedy hill-climber , i.e.
we tried to learn additional edges from the data, which are missing in the literature network.
this would be the intended and devised approach for using depns practice, if possible.
like before, we applied the nonparametric stratified bootstrap to derive edge probabilities.
as shown in figure  <dig> the initialization with known edges from the literature made all three spurious edges dispensable to explain our data.
two edges were inferred in addition to the initial network: pakt <dig> → perk1/ <dig> and cdk4: → cyclind <dig>  as mentioned before, pakt <dig> → perk1/ <dig> can be explained as an indirect activation of perk1/ <dig> by akt <dig> through igf1r.
furthermore, it has been observed that cdk <dig> and cyclind <dig> can form a complex phosphorylating prb  <cit> , which explains cdk <dig> → cyclind <dig>
in conclusion of our experiments we have shown that depns were able to reconstruct the erbb-regulated g1/s transition network in a sufficient quality, even if there were many latent nodes.
the approach presented in this paper differs significantly from the one in our earlier publication on this data set  <cit> , in which our aim was to explore biological implications of the literature network.
this was done by comparing outcomes of knock-down simulations to our experimental data.
specifically, in our previous paper no likelihood model of the data was employed and therefore no network inference could be made.
both of these complementary modeling approaches, i.e.
network reconstruction  and network based prediction of outcomes  <cit>  demonstrate a good consistency of the literature network with the experimental data we have.
related work
in this paper we introduced deterministic effects propagation networks as a special instance of bayesian networks for robustly estimating signaling networks from multiple and combinatorial intervention effects.
depns have certain similarities but also major differences to nested effects models   <cit> .
nems were especially designed to estimate signaling networks from high dimensional secondary intervention effects monitored on microarray data.
the idea is to reverse engineer upstream signaling cascades from the nested structure of downstream effects, which is observed from several single gene interventions.
the network is estimated between all perturbed genes.
although our data in this paper looks relatively different at a first glance , there are also some similarities: if a certain protein i is perturbed, we expect to see a perturbation of all downstream proteins as well.
that means the set of proteins being perturbed, if i is perturbed is a superset of the proteins being perturbed if any downstream protein of i is perturbed.
in fact with our approach we implicitly learn this nested effects structure from data and encode it into a network.
it should be mentioned that both, nems and depns, can be interpreted within the bayesian networks framework in a similar way  <cit> .
an additional feature of depns compared to current nem implementations is the possibility for latent network nodes.
the advantage of this property is a more comprehensive view of the network of interest.
however, it should be mentioned that it is only possible to estimate edges leading from latent nodes to measured nodes from data.
another major difference to nems is that depns reconstruct the network between affected nodes from combinatorial perturbations of the same nodes.
therefore, depns and nems may be seen as orthogonal approaches.
depns are quite different from other bayesian network approaches for learning networks from intervention effects  <cit> : depns do not assume a direct conditional dependence between measurement distributions.
instead depns rely on a deterministic effects propagation scheme, which makes all measurement distributions statistically independent.
therefore, depns can learn graphs structures containing loops, whereas usual bayesian networks are restricted to directed acyclic graphs.
this is also an important difference to  <cit> , which both employ ideal interventions  <cit> .
these approaches assume to have few perturbations and a large amount of measurements  without interventions.
in contrast, depns make the opposite assumption: perturbations for almost all proteins in the network, but only few measurements  without interventions.
in contrast to mechanistic, quantitative approaches  <cit> , depns are probabilistic network models.
depns do not aim to model real protein concentration changes over time , but purely make probabilistic inference on the network structure given experimental data.
depns differ significantly from reconstruction methods that overlay experimental data with known pathway information .
in contrast depns make a statistical inference of network structures from experimental measurements using a likelihood model.
CONCLUSIONS
studies of intervention effects play an increasing role in the understanding of complex biological networks.
we introduced deterministic effects propagation networks  as a bayesian network computational approach to infer non-transcriptional signaling networks from multiple and combinatorial low dimensional intervention effects.
in this paper these were quantitative proteomics data measured via reverse phase protein arrays.
our method hence aims at a system's understanding on the protein level.
nonetheless, depns are in principle not limited to proteomics data and could be applied to other data showing similar characteristics.
unlike nested effects models  <cit> , our method is especially designed for cases, where the number of perturbations exceeds the number of proteins in the network.
our method allows for the existence of latent nodes and of missing values.
the reconstruction accuracy of depns degrades with an increasing number of latent nodes within an expected range.
this phenomenon is up to a certain degree unavoidable and common to all statistical approaches, since latent nodes imply a complete lack of data for one or several variables, which need to be deduced from available information.
apart from that our method was shown to be highly robust in our simulation studies.
this specifically includes a very robust behavior against missing values.
we applied depns to reconstruct the erbb signaling network in trastuzumab resistant breast cancer cells from  <dig>  partially combinatorial, sirna interventions.
the data were measured via reverse phase protein arrays.
although the data contained a high fraction of nodes without measurements, our network reconstruction was highly specific and showed only three spurious interactions.
these could be made dispensable by starting the network reconstruction from current literature knowledge.
altogether our approach was able to reliable uncover a significant amount of the currently known erbb-regalated human g1/s transition network, even with relatively little data.
in conclusion of our work we believe that depns offer a robust and reliable approach to reverse engineer non-transcriptional signaling cascades from multiple low dimensional intervention effects.
compared to other bayesian networks, which model directly conditional dependencies between measurements, depns showed a much higher network reconstruction accuracy in our simulation studies.
at the same time depns have a simpler inference mechanism, which allows for learning cyclic network structures in a straight forward way.
a limitation of depns lies in the supposed deterministic way of effects propagation.
moreover, there is no inference made on the combinatorial functions carried out by the individual network nodes.
furthermore, depns currently do not offer a model of the time dynamics .
while a model of the time behavior made no sense for our specific data, it is worthwhile to investigate dynamic depns for different data types in the future.
clearly, there are lots of other possibilities to improve the current method in the future.
besides inference on combinatorial functions of network nodes, this includes a more efficient way to deal with missing values and latent nodes, e.g.
via structural em  <cit> , inference of up- and downregulated edges and improved network structure search strategies  <cit> .
these issues should be investigated in the future.
the complete method as well as the data set is available as part of the r package "nem"  <cit>  as a supplement to this paper .

 reverse engineering of biological networks is a key for the understanding of biological systems. the exact knowledge of interdependencies between proteins in the living cell is crucial for the identification of drug targets for various diseases. however, due to the complexity of the system a complete picture with detailed knowledge of the behavior about the individual proteins is still in the far future. nonetheless, the advent of gene perturbation techniques, like rna interference, has enabled the possibility to study cellular systems systematically under varying conditions, hence opening new perspectives for network reconstruction methods. a number of approaches have been proposed in the literature for estimating networks from perturbation effects. many of these approaches aim at reconstructing a network from directly observable effects. for example, rung et al. builds a directed disruption graph by drawing an edge, if gene i results in a significant expression change at gene j. wagner uses such disruption networks as a starting point for a further graph-theoretic method, which removes indirect effects, hence making the network more parsimonious. tresch at el. extend this approach by additionally making use of p-values and fold-change directions to make the network more consistent with the observed biological effects. also bayesian networks have been used to model the statistical dependency between perturbation experiments. for this purpose pearl proposes an idealized model of interventions. he assumes that once a network node is manipulated, the influence of all parent nodes is eliminated and the local probability distribution of the node becomes a point mass at the target state. besides for bayesian networks, ideal interventions have also been applied for factor graphs and dependency networks. epistasis analysis offers a possibility for learning from indirect downstream effects. for example, driessche et al. use expression profiles from single and double knockdowns to partly reconstruct a developmental pathway in d. discoideum via a simple cluster analysis. also fully quantitative models using differential equation systems have been suggested. for example, nelander et al. propose a model for predicting combinatorial drug treatment effects in cancer cells. recently, nested effects models have been proposed as a method, which is specifically designed to learn the signaling flow between perturbed genes from indirect, high-dimensional effects, typically monitored via dna microarrays. nems use a probabilistic framework to compare a given network hypothesis with the observed nested structure of downstream effects. perturbing one gene may have an influence on a number of downstream genes, while perturbing others affects a subset of those. moreover, several of these subsets could be disjoint, i.e. the knockdown of gene i shows effects, which mostly do not overlap with the effects seen at the knockdown of gene j. nems have been applied successfully to data on immune response in drosophila melanogaster, to the transcription factor network in saccharomices cerevisiae, to the er-α pathway in human breast cancer cells, and to the synthetic lethality interactions network in saccharomicies cerevisiae. the work presented in this paper is designed for a different scenario: we would like to reverse engineer a protein signaling network based on experimentally measured effects on protein expression and activation level after multiple interventions. these interventions may also be combinatorial, i.e. there is more than one knock-down at a time. importantly, the set of all perturbations should cover a fraction as large as possible of the network proteins. effects of all interventions on the network proteins are quantified directly on protein expression and activation level, for instance via reverse phase protein arrays. here, we propose a probabilistic approach called deterministic effects propagation networks, which can estimate the most likely signaling network based on these data. depns are a special case of bayesian networks, which employ a mixture of purely deterministic and gaussian variables. while depns and nems have a similar effects propagation scheme, depns differ from nems with respect to the following aspects: • in nems each node corresponds to one perturbation experiment. in depns each node corresponds to one single protein, potentially influenced by one or several perturbations. • nems only work, if the number of perturbations is much smaller than the number of monitored downstream effects. in depns the opposite is true: depns assume that at least each network node has been perturbed once. more perturbations are beneficial. • in nems indirect, high dimensional downstream intervention effects are monitored. in depns intervention effects are monitored on all other network proteins. • depns offer the possibility to include combinatorial interventions in addition to single interventions in a natural way. • depns can use latent network nodes. this has not been introduced for nems so far. • depns offer a mechanism for missing value imputation. this has not been introduced for nems so far as well. in contrast to other bayesian networks, depns allow for directed graph structures with loops, because they employ a deterministic effects propagation scheme. our approach was first validated using an extensive simulation study. afterwards it was employed to reverse engineer the erbb receptor regulated g1/s transition network in hcc human breast cancer cells, which show de novo resistance against trastuzumab. in this dataset knock-downs of proteins utilizing rnai were performed and effects were measured via rppas.  deterministic effects propagation networks our basic assumption is that we measure an unknown signaling network of proteins by performing multiple interventions on these proteins. these perturbations may affect single network proteins at a time or may also be combinatorial, i.e. they affect several network proteins at one time. for each network protein in the network we monitor the effect of all interventions carried out on proteins in the network. note, that each specific intervention not only influences direct targets, but may also cause effects on downstream network proteins. we explicitly assume that in one experiment all network proteins are unperturbed. we employ special instances of bayesian networks, which we call deterministic effects propagation networks, to infer the most likely protein network given measurements from multiple interventions under different experimental conditions or time points. our implementation of depns can deal with latent network nodes and with missing data. the details of our method are described in the section material and methods. briefly, the idea is that we have an unknown network graph, where each node can have two states: perturbed and unperturbed. in principle, this could also be extended to three states, but in this paper we only deal with the simpler binary case. furthermore, attached to each node we have experimental measurements, which for each time point are assumed to come from two gaussian distributions, one for the perturbed and one for the unperturbed case. the parameters of these distributions are determined in the following way: for each perturbation experiment we compute the expected downstream effects. this is done by assuming node i to be perturbed whenever i itself or any of its parents are perturbed, i.e. perturbation effects are propagated in a deterministic fashion from parents to children. this implies that the network structures learned by depns are always transitively closed graphs: if there is a path from a to b then there is also an edge a → b. it should to mentioned that the perturbations state of a protein can never depend on itself through a loop in our model, i.e. when the state of a variable is obtained it is not updated further. however, this behavior could be changed in principle. since we have several perturbation experiments and also one control experiment, we can now determine for each protein measurement, whether the protein was perturbed or not. with this information the data for each protein gets divided into two distributions: one for the perturbed and one for the unperturbed case. the parameters of these distributions, which are supposed to be gaussian, are estimated either in a maximum likelihood or in a bayesian fashion. the details are described in the material and methods section. given the parameters we can then calculate the probability of observing the data for all proteins under a given network hypothesis. an edge a → b in this network denotes an influence of a on b. this may be interpreted as the direction of the signal flow. our model does not distinguish between signaling via a transcriptional regulation or via a protein phosphorylation. the type of the edge is dependent on the type of antibody used. if the number of nodes in the network is considerably small, we can enumerate all possible network hypothesis and take the one with the highest likelihood. however, this exhaustive search will be impossible as soon as the network gets larger. for instance, for nodes there are already network hypothesis to test. therefore, here we restrict ourselves to a greedy hill climbing approach, in which we begin with an initial network and then successively add the edge increasing the likelihood most. all possible alternatives are checked. we do not use any arbitrary starting point for the hill climber, because depns operate in the space of transitively closed graphs only. therefore, no backtracking/edge removal step can be implemented in an obvious way. hence, initial networks should be chosen in a sensible way. simulation studies graph sampling we tested the performance of our approach in a simulation study on artificial data, which were created as follows: subgraphs with n = and n = nodes were cut out randomly from randomly selected kegg signal transduction pathways. only gene-gene interactions were considered as edges in the graph. please note that our sampled subgraphs are not random graphs, but randomly selected existing subgraphs of kegg pathways. for n = % of the sampled subgraphs were acyclic. for n = these were %. data simulation single knock-downs of all individual network proteins as well as one experiment without interventions were simulated as follows: for each simulated experiment m replicates were drawn from normal distributions for each network node. the means of these normal distributions were sampled from n, if the node was expected to be perturbed according to the deterministic effects propagation scheme described above, and from n otherwise. the variances of the normal distributions were drawn from inv-χ the choice of these parameters corresponds to those in section "methods" and is in agreement with our experimental data. network reconstruction figures and visualize the performance of a network reconstruction with our approach for a varying number of replicates in terms of sensitivity and specificity. for each number of replicates networks with simulated perturbation data were generated as described before. as seen from the plots our method achieves a constantly high specificity of >90% with a sensitivity between 70% and 80%, i.e. the number of replicates has only a minor influence on the quality of the network reconstruction. the reason for this robust behavior is that due to the effects propagation scheme, we only operate in the space of transitively closed graphs. the restricted model class can be viewed as a kind of regularization helping to identify a well performing network structure. additionally, effects propagation increases the chance of getting enough data for both, the perturbed and the unperturbed case, for reliable parameter estimation. we compared our depn approach to a standard implementation of bayesian networks, which directly models the dependency of measurements: for this purpose we discretized our simulated data. this was done by setting if the data entry was closer to than to, and otherwise. the conditional distribution at each network node was therefore a binomial with parameter equal to the perturbation probability of the node. we employed the r package bnlearn to learn the bn network structure via a hill climbing procedure initialized with an empty network. the bayesian dirichlet score was used to evaluate network structures. please note, that bnstandard in contrast to depns can only learn network structures without cycles. however, the vast majority of all original kegg graphs were acyclic here. since our simulated data contains only single knock-outs, we also compared our depn approach to nested effects models: for this purpose we discretized our simulated data as before. the likelihood model introduced in was employed with type i and type ii error rates, which are expected by the known variances and means of the perturbed and unperturbed distributions. the same greedy hillclimber as for depns was employed for network structure search. figures and show a similar performance of bnstandard compared to depns in terms of specificity, but a significantly worse sensitivity. the performance on acyclic graphs was essentially the same as on all graphs. in terms of computation times bnstandard is clearly faster than the depn approach. future optimizations of our depns implementation may reduce this gap. nevertheless, depns are only suitable for small network reconstructions. this is, however, in the light of the large experimental effort for multiple interventions and determination of their respective effects on network proteins not really a practical limitation. in summary the comparison clearly demonstrates that a standard bayesian network approach, which models only conditional dependencies between measurements while ignoring the perturbation signal flow in the network, is not suitable for our kind of data. compared to nems, depns show a significantly higher specificity and a slightly reduced sensitivity. in terms of computation times, nems are clearly faster than depns. this is expected, because nems can be formulated in a much more efficient way. since nems only employ single gene perturbations, they allow to perform deterministic effects propagation in an implicit rather than explicit way and do not have to estimate distribution parameters like depns. apart from that nems and depns have a lot of similarities, but also major differences, which we highlight in section "related work". altogether, both comparisons demonstrate that depns seem to be well suited for the low dimensional knock-down effects data, which we are interested in here. missing values we next investigated the effect of an increasing number of missing values in a dataset with replicate measurements in networks with n = and n = nodes, i.e. × × = and × × = data points in total. the positions of missing values were chosen uniform randomly in the complete data matrix. as seen from figures our method behaves very robust against replacing missing values with estimated posterior distribution mode values, even with out of and out of missing values in total. again, this robust behavior can be attributed to the same reasons as described before: a restricted model class together with the deterministic effects propagation scheme. latent nodes we also monitored the effect of having latent nodes in the network. again, data with replicates was simulated in networks with n = and n = nodes. the results shown in figures reveal a constantly high specificity of >90%, whereas the sensitivity drops with latent node from ~80% to ~60% with n = and ~80% to ~70% with n = respectively. with latent nodes we still achieve a sensitivity of ~40% with n = and > 60% with n = with latent out of nodes we have ~20% sensitivity, and with latent out of nodes we have ~40% sensitivity. we have to remind here that only outgoing edges from latent nodes can be estimated, hence the result shows an expected behavior. note that on average with each latent node the number of unobservable edges should increase by roughly 1/ = ~17% with n = and by 10% with n =. we believe that a positive point of our method is that still the specificity is rather high, meaning we do not get too many false positive edges. reconstruction of the erbb receptor-regulated g1/s transition network rnai and reverse phase protein arrays sixteen proteins playing a role in the g1/s transition of human cells were taken. among these proteins there were the receptors erbb erbb and erbb several signaling intermediates and cell cycle proteins. finally, we included the retinoblastoma protein prb which regulates the g1/s transition. interactions between these proteins known from the literature are shown in figure and additional file sixteen knock-downs with chemically synthesized sirnas for proteins in our network and one experiment with mock transfected cells were conducted. after each sirna transfection cells were stimulated with egf for hours. we measured the protein expression of signaling intermediates for each knockdown sample using reverse phase protein arrays before and after egf stimulation with technical and biological replicates. these data were integrated using quantile normalization. further details of the experimental setup and a list of performed interventions can be found in the material and methods part of this paper. additional information can be found in our earlier publication. a heatmap of the data matrix is available in the supplement of this paper. network reconstruction with depns we used depns to learn a network between the sixteen proteins purely from data. since for out of the proteins in our network no measurements were available due to the lack of antibodies, they were marked as latent nodes. furthermore, for perk1/ we had only instead of biological replicates for cyclin e and cyclin d knock-downs. since our method in the current implementation requires the same amount of data for each experiment, measurements for the third biological replicate were treated as missing values and imputed during network learning using the em algorithm approach described in the materials and methods part of this paper. network structure learning was performed via greedy hill climbing. in order to retrieve a measure of confidence for the individual edges learned by our method, we performed a nonparametric stratified bootstrap. that means we sampled times with replacement from each replicate group within each condition. on each bootstrap sample we estimated a network using depns. after the bootstrap we counted the relative frequency of each edge in all depn reconstructions. the relative frequency is an approximate probability for each edge. only edges with a probability > were considered as being reliable. the result shown in figure reveals that our methods correctly recovered 20/ interactions from the literature network with probability >, although we had latent nodes. here we also counted edges as correct, if they correspond to indirect interactions in the original network. this was valid, since depns learn transitively closed graphs. the edge erbb → erbb could be explained by the fact that erbb and erbb together form a heterodimer activating the mapk and pi3k pathway. likewise, the edge cdk → cdk was there, because both proteins form a complex. the complex formation of cdk and cdk could also explain the edge cdk → p. it was found that akt increases the expression of human igf1r, which explained the edge pakt → perk1/ as an indirect activation of perk1/ by akt through igf1r. in conclusion there were only three edges, namely cdk → erbb pakt → erbb and myc → erbb for which we could not find any explanation in the literature. we additionally performed a network inference starting from the literature network as an initial network for the greedy hill-climber, i.e. we tried to learn additional edges from the data, which are missing in the literature network. this would be the intended and devised approach for using depns practice, if possible. like before, we applied the nonparametric stratified bootstrap to derive edge probabilities. as shown in figure the initialization with known edges from the literature made all three spurious edges dispensable to explain our data. two edges were inferred in addition to the initial network: pakt → perk1/ and cdk4: → cyclind as mentioned before, pakt → perk1/ can be explained as an indirect activation of perk1/ by akt through igf1r. furthermore, it has been observed that cdk and cyclind can form a complex phosphorylating prb, which explains cdk → cyclind in conclusion of our experiments we have shown that depns were able to reconstruct the erbb-regulated g1/s transition network in a sufficient quality, even if there were many latent nodes. the approach presented in this paper differs significantly from the one in our earlier publication on this data set, in which our aim was to explore biological implications of the literature network. this was done by comparing outcomes of knock-down simulations to our experimental data. specifically, in our previous paper no likelihood model of the data was employed and therefore no network inference could be made. both of these complementary modeling approaches, i.e. network reconstruction and network based prediction of outcomes demonstrate a good consistency of the literature network with the experimental data we have. related work in this paper we introduced deterministic effects propagation networks as a special instance of bayesian networks for robustly estimating signaling networks from multiple and combinatorial intervention effects. depns have certain similarities but also major differences to nested effects models. nems were especially designed to estimate signaling networks from high dimensional secondary intervention effects monitored on microarray data. the idea is to reverse engineer upstream signaling cascades from the nested structure of downstream effects, which is observed from several single gene interventions. the network is estimated between all perturbed genes. although our data in this paper looks relatively different at a first glance, there are also some similarities: if a certain protein i is perturbed, we expect to see a perturbation of all downstream proteins as well. that means the set of proteins being perturbed, if i is perturbed is a superset of the proteins being perturbed if any downstream protein of i is perturbed. in fact with our approach we implicitly learn this nested effects structure from data and encode it into a network. it should be mentioned that both, nems and depns, can be interpreted within the bayesian networks framework in a similar way. an additional feature of depns compared to current nem implementations is the possibility for latent network nodes. the advantage of this property is a more comprehensive view of the network of interest. however, it should be mentioned that it is only possible to estimate edges leading from latent nodes to measured nodes from data. another major difference to nems is that depns reconstruct the network between affected nodes from combinatorial perturbations of the same nodes. therefore, depns and nems may be seen as orthogonal approaches. depns are quite different from other bayesian network approaches for learning networks from intervention effects: depns do not assume a direct conditional dependence between measurement distributions. instead depns rely on a deterministic effects propagation scheme, which makes all measurement distributions statistically independent. therefore, depns can learn graphs structures containing loops, whereas usual bayesian networks are restricted to directed acyclic graphs. this is also an important difference to, which both employ ideal interventions. these approaches assume to have few perturbations and a large amount of measurements without interventions. in contrast, depns make the opposite assumption: perturbations for almost all proteins in the network, but only few measurements without interventions. in contrast to mechanistic, quantitative approaches, depns are probabilistic network models. depns do not aim to model real protein concentration changes over time, but purely make probabilistic inference on the network structure given experimental data. depns differ significantly from reconstruction methods that overlay experimental data with known pathway information. in contrast depns make a statistical inference of network structures from experimental measurements using a likelihood model.  studies of intervention effects play an increasing role in the understanding of complex biological networks. we introduced deterministic effects propagation networks as a bayesian network computational approach to infer non-transcriptional signaling networks from multiple and combinatorial low dimensional intervention effects. in this paper these were quantitative proteomics data measured via reverse phase protein arrays. our method hence aims at a system's understanding on the protein level. nonetheless, depns are in principle not limited to proteomics data and could be applied to other data showing similar characteristics. unlike nested effects models, our method is especially designed for cases, where the number of perturbations exceeds the number of proteins in the network. our method allows for the existence of latent nodes and of missing values. the reconstruction accuracy of depns degrades with an increasing number of latent nodes within an expected range. this phenomenon is up to a certain degree unavoidable and common to all statistical approaches, since latent nodes imply a complete lack of data for one or several variables, which need to be deduced from available information. apart from that our method was shown to be highly robust in our simulation studies. this specifically includes a very robust behavior against missing values. we applied depns to reconstruct the erbb signaling network in trastuzumab resistant breast cancer cells from partially combinatorial, sirna interventions. the data were measured via reverse phase protein arrays. although the data contained a high fraction of nodes without measurements, our network reconstruction was highly specific and showed only three spurious interactions. these could be made dispensable by starting the network reconstruction from current literature knowledge. altogether our approach was able to reliable uncover a significant amount of the currently known erbb-regalated human g1/s transition network, even with relatively little data. in conclusion of our work we believe that depns offer a robust and reliable approach to reverse engineer non-transcriptional signaling cascades from multiple low dimensional intervention effects. compared to other bayesian networks, which model directly conditional dependencies between measurements, depns showed a much higher network reconstruction accuracy in our simulation studies. at the same time depns have a simpler inference mechanism, which allows for learning cyclic network structures in a straight forward way. a limitation of depns lies in the supposed deterministic way of effects propagation. moreover, there is no inference made on the combinatorial functions carried out by the individual network nodes. furthermore, depns currently do not offer a model of the time dynamics. while a model of the time behavior made no sense for our specific data, it is worthwhile to investigate dynamic depns for different data types in the future. clearly, there are lots of other possibilities to improve the current method in the future. besides inference on combinatorial functions of network nodes, this includes a more efficient way to deal with missing values and latent nodes, e.g. via structural em, inference of up- and downregulated edges and improved network structure search strategies. these issues should be investigated in the future. the complete method as well as the data set is available as part of the r package "nem" as a supplement to this paper. 
depns allow to reconstruct protein networks based on combinatorial intervention effects, which are monitored via changes of the protein expression or activation over one or a few time points. in this paper, we propose deterministic effects propagation networks as a special bayesian network approach to reverse engineer signaling networks from a combination of protein expression and perturbation data. depns offer a robust, efficient and simple approach to infer protein signaling networks from multiple interventions. we applied depns to reconstruct the erbb signaling network in de novo trastuzumab resistant human breast cancer cells, where protein expression was monitored on reverse phase protein arrays after knockdown of network proteins using rnai. in combination with mrna or protein expression data this allows to gain insights into the behavior of complex biological systems. our implementation of depns allows for latent network nodes and has a built in mechanism to impute missing data. the method as well as the data have been made part of the latest version of the r package "nem" available as a supplement to this paper and via the bioconductor repository.
 modern gene perturbation techniques, like rna interference, enable us to study effects of targeted interventions in cells efficiently. in combination with mrna or protein expression data this allows to gain insights into the behavior of complex biological systems.  in this paper, we propose deterministic effects propagation networks as a special bayesian network approach to reverse engineer signaling networks from a combination of protein expression and perturbation data. depns allow to reconstruct protein networks based on combinatorial intervention effects, which are monitored via changes of the protein expression or activation over one or a few time points. our implementation of depns allows for latent network nodes and has a built in mechanism to impute missing data. the robustness of our approach was tested on simulated data. we applied depns to reconstruct the erbb signaling network in de novo trastuzumab resistant human breast cancer cells, where protein expression was monitored on reverse phase protein arrays after knockdown of network proteins using rnai.  depns offer a robust, efficient and simple approach to infer protein signaling networks from multiple interventions. the method as well as the data have been made part of the latest version of the r package "nem" available as a supplement to this paper and via the bioconductor repository. 