BACKGROUND
proteins’ rates of evolution vary across orders of magnitude .
understanding the factors underlying this enormous variability is a central problem in evolutionary biology.
the neutral theory of molecular evolution predicts higher selective pressures  acting on genes encoding proteins that perform the most important functions for organisms’ biological fitness
<cit> ; however, the relative contributions of genes to fitness remain elusive.
over the last decade, the availability of both genomic and functional data has allowed us to identify some of the factors that correlate with genes’ rates of evolution, among which patterns and levels of gene expression appear to be the most important .
however, this gene-centered approach has been able to explain only a relatively small fraction of the variability of rates of evolution .
genes and proteins do not act in isolation, but usually operate as components of complex networks of interacting molecules.
therefore, taking into account the interaction networks in which a gene participates may provide key insight into the evolutionary forces acting on it.
furthermore, understanding the function and evolution of molecular networks can aid applications such as metabolic engineering and drug discovery and design .
molecular networks in general can be modeled as undirected graphs, whose nodes and edges represent proteins and interactions, respectively.
over the last years, the accumulation of interactomic data has allowed us to cast a first glance at the impact of genes’ position in the network on their patterns of molecular evolution .
for instance, genes encoding the most connected proteins in metabolic and protein–protein interaction networks are more selectively constrained than those acting at the network periphery
<cit> , and interacting genes tend to show correlated evolutionary histories, e.g.
evolving at similar rates
<cit> .
these observations indicate that the structure of molecular networks impose constraints on the evolution of their components.
in addition to the centrality of proteins in functional networks, their hierarchical position in molecular pathways has also been put forward as a determinant of their rates of evolution.
remarkably, in a number of biosynthetic pathways, the strength of purifying selection has been found to correlate with the position of genes along the upstream/downstream pathway axis, with upstream genes being generally the most constrained
<cit> .
this distribution of the levels of selective constraint has been attributed to upstream genes being involved in the biosynthesis of a broader range of biochemical compounds than downstream genes, and hence being more pleiotropic
<cit> , and/or to upstream genes exerting a higher influence over the flux of metabolites along pathways
<cit> .
consistent observations have been made in transcriptional regulatory networks: upstream genes in these networks are more likely to be essential
<cit> , and simulation analysis has shown that genes exerting a higher degree of control over other genes, and/or those that are less regulated by other genes, are more strongly affected by selection
<cit> .
less well understood is how selective pressures distribute along signal transduction pathways.
these pathways play a key role in modulating cell function in response to extracellular and intracellular stimuli.
they are characterized by the presence of a receptor that, in response to signals, is able to trigger a cascade of events that ultimately modulate the final effectors, which are responsible for mediating the biological responses of the pathway.
the signal is transmitted from the upstream to the downstream part of the pathway by means of a series of molecular interactions .
these interactions have a directed nature, i.e., each interaction involves an upstream modulator  and a downstream modulated element .
therefore, signaling pathways and networks can be represented as directed graphs, in which nodes and arcs represent molecules and directed interactions, respectively.
despite the availability of some directed interaction networks , their directed nature has often been overlooked in evolutionary analyses, being often treated as undirected networks.
so far, the distribution of levels of selective constraint across the upstream/downstream axis of signaling pathways has been studied in only a few pathways, with contrasting results.
in the drosophila ras pathway, upstream genes appear to be more selectively constrained than those acting at the downstream part
<cit> .
this is to be expected under a model in which upstream genes, owing to their higher hierarchical position, have a higher impact on the function of the pathway.
indeed, a single activated upstream protein has the potential to activate several downstream target molecules ; therefore, higher levels of selective constraint might be expected at these genes.
however, the opposite pattern  has been described in the drosophila, caenorhabditis and vertebrate insulin/tor pathway
<cit> , and in the yeast hog pathway
<cit> .
in both the yeast hog pathway and the caenorhabditis insulin/tor pathway, the tendency appears to vanish when the effect of confounding factors that correlate with both network position and levels of selective constraint, such as expression level or codon bias, are removed
<cit> , whereas in the drosophila and vertebrate insulin/tor pathway the polarity in the levels of selective constraint is independent of these and other factors
<cit> .
therefore, the effect of upstream/downstream position of genes in signal transduction pathways on their rates of evolution, and the relative contribution of potential confounding factors, remain elusive.
here, i combine comparative genomics and manually-curated directed protein interaction data to determine whether purifying selection differentially acts on the upstream and downstream genes of the human signal transduction network.
for that purpose, new methods that explicitly take into account the directed character of the network are introduced.
surprisingly, no evidence was found supporting that upstream genes  are subject to increased levels of selective constraint.
on the contrary, purifying selection was found to act more severely on genes acting at the downstream part of the network, which seems to be due to downstream genes being more highly and broadly expressed, having certain functions and, in particular, encoding proteins that are more highly connected in the protein–protein interaction network.
when the effect of these confounding factors was discounted, upstream and downstream genes were found to evolve at similar rates.
these results clearly indicate that the hierarchical position of genes in the signal transduction network has, in general, no direct effect on their rates of evolution, suggesting that genes occupying more hierarchical positions in the network are not more relevant for the function of the network and the fitness of the organism.
the trend found in the overall signaling network is exemplified by analysis of the distribution of purifying selection along the mammalian ras signaling pathway, showing that, contrary to previously suggested in drosophila, upstream and downstream genes evolve at similar rates.
these observations imply that natural selection differently distributes across signal transduction networks and across biosynthetic and transcriptional regulatory networks.
this might reflect fundamental differences in the function and organization of the signal transduction network and metabolic and transcriptional regulatory networks.
RESULTS
unexpectedly, genes acting at the upstream part of the human signal transduction network exhibit higher rates of evolution
the human signal transduction network gleaned by cui et al.
<cit>  was used.
in this network, directed edges  between pairs of proteins represent regulatory interactions  — e.g., the directed edge “a → b” denotes that protein “a” activates or inhibits protein “b”.
the dataset was assembled and manually curated by cui et al.
by merging the interactions compiled in different databases, followed by validation of each interaction by accurate revision of the literature.
undirected edges, which represent generic  physical interactions, were not considered in the current analysis, since they contain no information about the hierarchical  position of genes in the network.
after filtering the dataset , it comprised  <dig> proteins connected by  <dig> non-redundant directed interactions .
this dataset included  <dig> ligands,  <dig> receptors,  <dig> adapters,  <dig> kinases,  <dig> phosphatases,  <dig> transcription factors,  <dig> structural proteins, and  <dig> ion transmembrane transporters.
most edges represent direct physical protein–protein interactions, with the only exception of a small fraction of “transcription factor → target gene” edges.
for each gene in the network, the 1: <dig> mouse ortholog was identified and the impact of natural selection was characterized from the nonsynonymous to synonymous divergence ratio .
values of ω equal to  <dig> are expected for genes evolving neutrally, whereas ω <  <dig> is indicative of the action of purifying selection preserving the sequence of the encoded proteins, and ω >  <dig> in a number of codons is indicative of the action of positive selection  driving the fixation of nonsynonymous substitutions.
the estimated ω values range from  <dig>  to  <dig>  , indicating that purifying selection has played an important role in the evolution of these genes.
i evaluated whether ω values depend on the hierarchical position that genes occupy in the network.
for that purpose, three different approaches were adopted to compare the evolutionary rates of upstream and downstream network genes.
first, the ω value of each gene was compared with the central value of the ω values of all its direct downstream targets .
out of the  <dig> genes in the network,  <dig> encode proteins that have at least one direct downstream target in the network .
of those,  <dig> genes exhibit a higher ω value than their targets, whereas  <dig> show lower ω values than their targets.
this asymmetry is not compatible with a random distribution of ω values across upstream and downstream genes .
therefore, surprisingly, upstream genes tend to show higher ω values than the genes that they modulate.
in order to establish whether the evolution of nonsynonymous or synonymous sites is responsible for this trend, this analysis was conducted for dn and ds separately.
clearly significant differences were observed for dn , but only marginally significant differences were observed for ds , indicating that the evolutionary patterns of nonsynonymous sites are the main responsible for the trend observed in ω.
*, p< <dig> ; **, p< <dig> , ***, p< <dig> .
anumber of genes with at least one direct downstream target .
only genes with available information for the parameter of interest were considered.
bnumber of genes with a higher value than the central value of its downstream targets.
cnumber of genes with a lower value than the central value of its downstream targets.
next, evolutionary rates of genes encoding proteins occupying extreme upstream and downstream positions in the network were compared.
a total of  <dig> genes were classified as upstream  and  <dig> were classified as downstream .
again, upstream genes were found to exhibit significantly higher ω values: upstream and downstream genes exhibit median ω values of  <dig>  and  <dig> , respectively .
significant differences were also observed in the dn , but not in the ds  values.
p-values were obtained from the mann–whitney u test.
*, p <  <dig> .
finally, i evaluated whether ω values correlate with parameters that reflect the hierarchical position of genes in the network.
both in-degree  and out-degree  negatively correlate with ω values, consistent with previous observations that highly connected genes are subject to stronger purifying selection
<cit> .
remarkably, the correlation is stronger for in-degree  than for out-degree , indicating that the ω value of a gene is more influenced by the number of upstream genes , than by the number of downstream targets .
similar results were obtained for dn, whereas ds does not correlate with in-degree or out-degree .
beyond considering the direct downstream targets of a gene, i also considered the effect of the total number of genes acting directly or indirectly downstream of a given gene  on its evolutionary rate.
this parameter positively correlates with ω , again indicating that proteins potentially affecting the activity of a higher number of other proteins  are encoded by genes that exhibit higher ω values.
again, significant results were obtained for dn , but not for ds .
*, p <  <dig> ; ***, p <  <dig> .
taken together, these results show that evolutionary rates exhibit a polarity across the human signaling network, with downstream genes progressively exhibiting lower ω and dn values.
this observation is in sharp contrast with patterns observed in metabolic and transcriptional regulatory pathways, in which upstream genes are generally the ones exhibiting lower rates of evolution
<cit> .
i next considered a number of potential factors that could account for the observed polarity.
the higher rates of evolution of upstream genes are the result of lower levels of purifying selection rather than a higher incidence of positive selection
the lower ω and dn values observed in the downstream network genes suggest a higher strength of purifying selection acting on these genes.
however, ω and dn values are also affected by the action of positive selection, and hence the observed trends could alternatively be the result of a putatively higher incidence of positive selection on upstream genes.
in order to discard this possibility, analyses were repeated after eliminating genes with the signature of positive selection.
when applied to pairs of human-mouse orthologs, the m <dig> vs. m <dig> test
<cit>  identified  <dig> network genes with potential signatures of positive selection .
excluding these genes yielded similar results.
first, after excluding these genes,  <dig> genes have at least one direct downstream target in the network, out of which  <dig> and  <dig> show higher and lower ω values than their targets, respectively .
significant results were also obtained for dn .
second, the  <dig> genes classified as upstream that display no signals of positive selection exhibit significantly higher ω , and dn values  than the  <dig> downstream genes without the signature of positive selection.
finally, ω and dn values positively correlate with h , and more strongly correlate with in-degree  than with out-degree .
because detection of positive selection is sensitive to the number of sequences used
<cit> , analyses were repeated using a total of  <dig> mammalian genomes .
a total of  <dig> network genes exhibit a 1: <dig> orthology across all genomes, out of which  <dig> exhibit potential signals of positive selection.
again, excluding these genes did not alter the observations that genes acting at the downstream part of the network exhibit significantly lower ω  and dn values .
therefore, the distribution of ω and dn values across the network may be the result of a polarity in the strength of purifying selection, with downstream genes being subject to higher levels of selective constraint, and not the result of a higher incidence of positive selection in the upstream part of the network.
subcellular location does not account for the higher rates of evolution of upstream network proteins
genes acting at the extracellular compartment have been reported to exhibit high rates of evolution
<cit> .
consistently, among genes in the used dataset, those acting at the extracellular space exhibit higher ω and dn values than the rest of the genes in the signal transduction network .
this could potentially explain the observed polarity in ω and dn values if proteins occupying highly hierarchical positions in the signaling network preferentially located to the extracellular space.
however, i tested this possibility by comparing the proportion of extracellular proteins among proteins acting at the upstream  and downstream  parts of the network, finding no significant differences .
i.e., proteins occupying extreme upstream positions in the network do not tend to act at the extracellular space, thus making it unlikely that the observed high ω values in genes encoding the upstream part of the network would be linked to subcellular location.
in order to further discard this possibility, analyses were repeated eliminating the  <dig> proteins targeted to the extracellular space, with similar results.
first, genes are significantly less selectively constrained than their downstream direct targets .
second, the  <dig> non-extracellular proteins classified as upstream exhibit higher ω and dn values than the  <dig> non-extracellular proteins acting at the downstream part of the network, although the differences are not statistically significant .
third, the ω and dn values correlate better with in-degree  than with out-degree .
although evolutionary rates do not significantly correlate with h, the direction of the correlations is positive and comparable in strength to that observed in the overall dataset .
these results indicate that the polarity observed in the levels of selective constraint is for the most part independent of the high evolutionary rates of proteins acting at the extracellular space.
using an expanded dataset , cui, purisima and wang
<cit>  found that signaling genes acting at the cell membrane, cytoplasm and nucleus exhibited different rates of evolution.
again, these differences raise the possibility that the observed polarity in the ω values across the network might be in part the result of the distribution of upstream and downstream proteins across the different subcellular compartments.
however, i found that  genes acting at the different cell compartments do not exhibit significantly different rates of evolution for the subset of genes included in the current analysis  , and that  the upstream and downstream parts of the network do not contain different proportions of proteins locating to the cell membrane , cytoplasm  or nucleus , making it unlikely that the polarity in the levels of selective constraint described here would obey to proteins’ subcellular location.
in order to further discard this possibility, analyses were conducted separately on subsets of genes acting at the different cell locations, with consistent results across all four compartments.
first, genes tend to evolve faster than their direct downstream targets acting at the same cell compartment, when the analyses are restricted to genes acting at the extracellular space , the cell membrane , the cytoplasm , or the nucleus .
second, ω values are higher for upstream than for downstream genes when genes acting at the extracellular space , the cell membrane , the cytoplasm  or the nucleus  are analyzed separately, although the differences are not statistically significant .
taken together, these results indicate that the polarity in the levels of selective constraint observed across the network is independent of the subcellular location of its components.
the lack of significance in some tests probably owes to the reduced statistical power resulting from partitioning the dataset.
gene functions upstream and downstream of the signal transduction network
genes with different functions are known to evolve at different rates .
hence, in order to gain further insight into the different selective pressures acting at the upstream and downstream parts of the human signaling network, it was tested whether the  <dig> genes acting at the upstream part performed different functions to the  <dig> genes acting at the downstream part.
for that purpose, the frequencies of the goslim goa terms of both groups was compared using the fatigo software
<cit> .
upstream genes exhibited a significant enrichment in the term “kinase activity” , whereas downstream genes were enriched in the category “ion transmembrane transporter activity” .
hence, the higher levels of selective constraint observed in the downstream part of the signaling network could potentially result from this unequal distribution of gene functions, if genes encoding kinases evolved particularly fast, and/or those encoding ion transmembrane transporters evolved particularly slow.
both gene groups were found to evolve slower than the rest of genes in the network , raising the possibility that the higher levels of selective constraint acting at the downstream part of the network would be due to the enrichment of this part of the network in ion transmembrane transporters.
in order to discard this possibility, analyses were repeated excluding the  <dig> network genes encoding ion transmembrane transporters, with similar results, although not statistically significant in all cases.
first, genes are significantly less selectively constrained than their downstream direct targets .
second, the  <dig> non-ion transporter proteins classified as upstream exhibit higher ω and dn values than the  <dig> non-ion transporter proteins acting at the downstream part of the network, although the differences are not statistically significant .
third, the ω and dn values correlate better with in-degree , than with out-degree , and they exhibit a positive, although non-significant, correlation with h .
when genes encoding kinases were also excluded from the analysis, the correlation between h and rates of evolution reached significance .
these results indicate that, although the different functions of upstream and downstream genes may partially contribute to the observed polarity in the levels of selective constraint, this factor does not completely account for the trend.
genes acting at the downstream part of the network exhibit higher expression levels and breadths, and encode more highly connected proteins
levels of selective constraint acting on a gene are affected by a number of factors, including gene expression level  and breadth
<cit> , codon bias
<cit> , connectivity of the encoded proteins
<cit> , and the number of paralogs
<cit> .
measures of these parameters were obtained from different sources .
because the dataset generated by cui et al.
<cit>  focuses on signaling molecules, proteins’ total connectivities were measured as the total number of physical interactors in the whole human protein–protein interaction network
<cit> .
all these factors significantly correlate with ω and dn in the dataset used .
therefore, a putative polarity in the distribution of these factors across the network could potentially account for the observed polarity in the levels of selective constraint.
remarkably, expression breadth and connectivity seem to be the best correlates of the rates of evolution, with both variables exhibiting a similar degree of correlation with ω and dn .
this result contrasts with previous observations that connectivity is a relatively weak predictor of rates of evolution , suggesting that this correlation is particularly strong for genes involved in signal transduction.
furthermore, the correlation between connectivity and rates of evolution remains significant when expression level and breadth, codon bias, and number of paralogs are simultaneously controlled for .
although genome-level analyses have revealed that the length of the encoded proteins also correlates with ω and dn <cit> , these correlations are not significant in the gene set used , and hence this parameter was not considered.
**, p< <dig> , ***, p< <dig> .
i evaluated whether factors correlating with ω and dn differ between upstream and downstream genes.
the paired test shows that genes tend to exhibit a lower expression breadth , and to encode proteins with a lower connectivity   than their direct downstream targets .
the test is only marginally significant for expression level .
as expected, both in-degree and out-degree positively correlate with connectivity.
the correlation is stronger for in-degree  than for out-degree , consistent with a higher connectivity in the downstream part of the network.
genes occupying extreme downstream positions in the pathway exhibit higher levels and breadths of expression, and encode more highly connected proteins , than upstream genes, although none of these differences is statistically significant .
the lack of significance in the mann–whitney tests might result from a reduced statistical power resulting from the fact that only a fraction of the genes were used in this analysis , whereas all genes were used in the other two analyses .
although none of the studied factors significantly depends on the upstream/downstream position of genes in all three analyses, these results provide evidence that expression breadth and connectivity, and perhaps expression level, exhibit a polarity across the human signaling network, with downstream genes being more highly and broadly expressed, and encoding more highly connected proteins than upstream genes.
highly and broadly expressed genes
<cit> , and those encoding highly connected proteins
<cit> , tend to be highly constrained.
combined, these trends raise the possibility that the polarity observed in the levels of selective constraint could be the result of the distribution of these factors across the network.
the lower evolutionary rates of downstream network genes are a byproduct of the distribution of protein connectivity and gene expression level and breadth across the network
i considered whether the observed polarity in the levels of selective constraint across the human signaling network was the result of the distribution of factors correlating with ω and dn.
for that purpose, this polarity was evaluated after discounting the putative effects of these factors.
linear regression was used to model the dependence of ω and dn from each factor separately, and the residuals of the model for each gene  were used in the paired and the mann–whitney tests.
the paired tests show that genes exhibit significantly higher ω and dn values than their downstream targets even if the effect of expression level, expression breadth, codon bias, or the number of paralogs is discounted, indicating that the difference in the intensity of purifying selection acting on upstream and downstream genes is independent of these factors.
however, when the effect of connectivity is discounted, the test yields no significant results .
genes occupying extreme upstream positions in the network are significantly less constrained than those acting at the downstream part, even when the effect of expression breadth, codon bias, or the number of paralogs is removed, but the difference is not significant when expression level or connectivity are included in the analysis .
given that the relationship between evolutionary rates and their correlates is not necessarily linear in all cases, a complementary nonparametric analysis was performed to assess the differences in ω and dn values between genes occupying extreme upstream and downstream positions.
partial correlation analysis was used to evaluate the association between upstream/downstream position  and rates of evolution while controlling for each of the controlling variables.
similar results were obtained: upstream and downstream genes exhibit significantly higher ω and dn values when controlling for codon bias or the number of paralogs, but not when controlling for expression level, expression breadth, or connectivity .
finally, partial correlation analysis was used to evaluate the association between h and ω and dn while controlling separately for each factor.
both correlations are significant when codon bias or the number of paralogs are used as controlling variables, but not when expression level, expression breadth, or connectivity are included in the analysis .
therefore, the differences in the ω and dn values of upstream and downstream network genes disappear when the effect of expression level, expression breadth, or connectivity, are removed.
this indicates that these factors account for the higher levels of selective constraint acting on the downstream genes of the human signal transduction network.
in particular, when connectivities are controlled for, ω and dn values are not significantly different for upstream and downstream genes, regardless of the method used , indicating that this factor may be the main responsible for the polarity in the levels of selective constraint.
when the effects of expression level or breadth are accounted for, the differences in ω and dn between upstream and downstream genes vanish for only one or two of the three methods used, indicating that these factors may contribute to the selective constraint polarity to a lesser extent.
a case study: natural selection upstream and downstream of the mammalian ras signaling pathway
so far, the drosophila ras pathway is the only signal transduction pathway for which higher levels of selective constraint have been described in genes acting at the upstream part
<cit> .
the ras pathway is activated by diverse extracellular stimuli, including growth factors, cytokines and hormones.
activation of receptor tyrosine kinases  triggers a cascade of signaling events  that ultimately promote cell proliferation, differentiation, migration and survival .
the architecture of this pathway has been characterized in detail in a number of organisms, including drosophila and mammals, revealing a highly conserved structure and function across metazoans.
riley et al.
<cit>  examined the patterns of molecular evolution of six of the genes involved in the drosophila ras pathway.
for that purpose, they surveyed the levels and patterns of polymorphism in drosophila melanogaster, and divergence with respect to the sister species d. simulans.
they found that three proteins acting at the upstream part of the pathway  were identical in both species, and exhibited no or just a few rare replacement polymorphisms within d. melanogaster.
on the contrary, the three downstream proteins included in their analysis  displayed a number of fixed differences between the two species, as well as some replacement polymorphisms in d. melanogaster.
consistent results were obtained from the ω values , pointing out to higher levels of purifying selection acting at the upstream part of the pathway.
a proposed explanation for this observation was that, owing to their highly hierarchical position, proteins acting at the upstream part of the pathway would act as key control points, greatly influencing the activity of downstream proteins, whereas downstream proteins would act as modulators of the output signal
<cit> .
the current wealth of genomic data provides an excellent opportunity to study the patterns of molecular evolution of all the core components of the ras pathway in other taxonomic groups.
hence, i considered whether levels of selective constraint acting on the mammalian ras pathway genes exhibit a similar distribution across the pathway to that observed in drosophila.
most of the afore-mentioned drosophila pathway genes have multiple orthologs in mammalian genomes known to participate in the ras pathway .
for each such mammalian gene, the ω value was estimated from comparison of the human-mouse pair of orthologs .
genes encoding the mammalian proteins sos <dig>  sos <dig>  erk <dig>  erk <dig>  and rsk1– <dig>  which are also key components of the pathway
<cit> , were also included in the current analysis .
contrary to previously suggested in drosophila <cit> , genes acting at the upstream part of the mammalian ras pathway  were not found to exhibit lower rates of evolution than those acting at the downstream part .
if anything, the opposite seemed to be the case , although the differences between both groups were not significant .
consistently, ω values were found to exhibit a negative, although non-significant, correlation with the position of genes along the upstream/downstream axis of the pathway .
for this analysis, pathway position was defined as the number of steps required to transduce the signal from the most upstream element , to each of the pathway components .
1: table s <dig> for a full list of the core pathway genes and their rates of evolution.
neither expression levels, nor expression breadths or connectivities significantly correlate with pathway position  or with ω values .
it should be noted, however, that microarray-based gene expression data, such as the one used here
<cit> , is subject to high levels of measurement noise.
more accurate measures of expression are available for some of the proteins involved in the ras pathway, and indeed it has been noticed that downstream proteins are present in higher concentrations in the cell , thus leaving open the possibility that the slightly higher levels of selective constraint observed in the downstream part of the pathway  may result from these genes being expressed at higher levels.
in any case, these results indicate that, contrary to previously suggested in drosophila <cit> , genes acting at the upstream part of the mammalian ras pathway are not subject to higher levels of selective constraint than those acting at the downstream part.
if anything, the opposite might be the case, although the trend is not statistically significant.
this observation is in agreement with the general trend described here for the whole human signal transduction network.
discussion
in summary, results presented here demonstrate that genes acting at the upstream part of the human signal transduction network clearly do not evolve under stronger selective pressures .
this finding seems striking, as these proteins, owing to their higher hierarchical position, may be capable of regulating the activity of a higher number of downstream proteins.
therefore, mutations in the genes encoding upstream network proteins would be expected to have more pleiotropic effects, and hence these genes would be expected to evolve under stronger purifying selection.
in contrast with these expectations, a number of analyses that capture different aspects of the hierarchical position of genes in the signaling network found no evidence for upstream genes being subject to higher levels of selective constraint.
on the contrary, upstream network genes are less selectively constrained than downstream genes, which appears to be explained by downstream genes performing certain functions, being more highly and broadly expressed and, in particular, encoding more highly connected proteins.
protein connectivity, gene expression breadth and, to a lesser extent, expression level, exhibit a polarity across the human signaling network, with downstream genes progressively encoding proteins that are more connected, and being expressed in a higher number of tissues and at higher levels .
the higher levels of expression of downstream genes might facilitate the amplification of signals through their progression throughout signaling cascades, which requires higher protein abundances for the downstream elements of the pathway .
additionally, the higher connectivity of downstream genes is consistent with previous observations that hubs  tend to be tightly regulated, exhibiting higher numbers of phosphorylation sites and higher mrna decay rates
<cit> .
the distribution of these factors across the network, combined with the fact that highly and broadly expressed genes
<cit> , and those encoding highly connected proteins in the protein–protein interaction network
<cit> , tend to evolve under stronger purifying selection, accounts for the higher levels of selective constraint observed in genes acting at the downstream part of the network.
in particular, when the effect of connectivity is factored out, the differences in ω and dn values between upstream and downstream genes are no longer significant, regardless of the method used , pointing to this factor as the main responsible for the higher levels of selective constraint acting on the downstream network genes.
when the effects of expression level and breadth are removed, the differences in the levels of selective constraint observed between upstream and downstream genes also vanish for some of the tests , suggesting that these factors also contribute to the tendency.
likewise, protein functions may also play a role, as the trend also vanishes for some tests when ion transmembrane transporters  are removed from the analyses.
in any case, taken together these observations indicate that the upstream/downstream position of proteins in the human signal transduction network does not have a direct effect on their rates of evolution.
the global trends observed in the entire human signal transduction network  mirror those previously observed in some particular pathways.
in both the caenorhabditis insulin/tor signaling pathway and the yeast hog signaling pathway, downstream genes evolve under stronger levels of selective constraint, probably as a result of their higher levels of expression and codon bias , respectively
<cit> .
it should be noted, however, that the general trends observed across the whole network do not necessarily imply that each particular signaling pathway must exhibit the same trend.
for instance, in the drosophila and vertebrate insulin/tor pathways, downstream genes evolve under stronger purifying selection, even when the putative effect of several factors that correlate with levels of selective constraint is discounted
<cit> .
another example that does not seem to conform to the general trend described here is the drosophila ras pathway
<cit> , whose upstream genes appear to be the most constrained, despite these genes being the ones expressed at lower levels
<cit>  .
therefore, despite the overall trends described here, the particular structure and function of each individual pathway may result in different distributions of the impact of evolutionary forces.
for instance, pathways whose upstream genes, for some reason, present higher connectivities, and/or levels or breadths of expression, would probably be more selective constrained in the upstream part.
it should also be noted that factors other than the ones considered here may also affect the distribution of levels of selective constraint across functional pathways and networks.
in particular, enzymes whose kinetic characteristics exert a high degree of influence over the entire network  are expected to be more selectively constrained than those exerting little influence on the network dynamics
<cit> .
the degree to which a particular enzyme can influence the behavior of the entire network is, however, difficult to determine from currently available data.
so far, the drosophila ras pathway is, to the best of the author’s knowledge, the only signal transduction pathway whose upstream genes appear to be the ones evolving under higher levels of purifying selection
<cit> .
however, due to the limited availability of genomic data at that time, this observation was based on the analysis of a small set of genes and the comparison of two closely related species .
results reported here for all core components of the mammalian ras pathway, and based on human-mouse comparisons , do not point in the same direction: upstream genes do not exhibit lower rates of evolution than those acting at the downstream part of the pathway.
if anything, the opposite might be the case .
in any case, this observation lends further support to the finding that genes occupying highly hierarchical positions in the human signal transduction network are not subject to higher levels of selective constraint than those acting at the downstream part.
the observation that genes acting at the upstream part of the human signal transduction network are not subject to higher levels of selective constraint sharply contrasts with the trends generally observed in biosynthetic pathways.
these pathways consist of a series of enzymes that share metabolites, with the product of a given enzyme being usually the substrate of the following enzyme  in the pathway.
in such pathways, metabolites can often take different alternative routes, which end up in the biosynthesis of different end products.
in a growing number of biosynthetic pathways, including the plant anthocyanin
<cit> , isoprene
<cit> , terpenoid
<cit> , and carotenoid pathways
<cit> , the bombyx melanin pathway
<cit> , and the mammalian dopamine pathway
<cit> , levels of selective constraint correlate with the position of genes along the upstream/downstream pathway axis, with upstream genes being always the most constrained.
two plausible scenarios might explain such a distribution of selective constraints in biosynthetic pathways.
first, as a result of their branching topology, these pathways often exhibit a hierarchical structure, with upstream genes being necessary for the biosynthesis of a higher number of end products than downstream genes; therefore, mutations in genes acting at the upstream part may have higher pleiotropic effects
<cit> .
second, upstream enzymes in biosynthetic pathways probably often exert a higher control over the flux of metabolites throughout the pathway
<cit> ; thus, the overall pathway function  may more strongly depend on the function of these enzymes.
both models predict a higher relevance for pathway function of upstream genes, and hence a higher strength of natural selection acting on these genes.
similar trends have been observed in transcriptional regulatory networks.
first, genes occupying higher hierarchical positions in these networks are likely to be essential, reflecting their key importance for organism’s fitness
<cit> .
second, simulation analyses of transcriptional regulatory networks have shown that genes exerting a high degree of control over other genes, or those that are less regulated by other genes, are expected to be more strongly affected by selection
<cit> .
therefore, results presented here imply that the distribution of evolutionary forces across the upstream and downstream parts of the signal transduction network remarkably differs from that observed in biosynthetic and transcriptional regulatory pathways.
in biosynthetic and transcriptional regulatory pathways, upstream genes evolve under stronger levels of selective constraint, whereas this is not the case for the human signal transduction network.
this suggests that genes occupying higher hierarchical positions are not more relevant for the function of the signaling network, or to the organism’s fitness, than downstream genes.
these different patterns of evolution may reflect fundamental differences in the function and organization of the different kinds of pathways.
as signaling networks become better understood, it becomes more apparent that signal transduction does not take place in the context of independent linear pathways in which each protein is exclusively controlled by the preceding protein in the cascade.
this classical model is being progressively substituted by a much more integrative model of signal transduction, in which a considerable amount of cross-talk takes place between pathways.
as a result, pathways would cooperate to integrate signals, and to make decisions about the responses in a decentralized manner, thus conferring robustness to the system
<cit> .
this kind of organization may allow the network to accommodate modifications in its upstream genes.
for instance, even if an upstream gene undergoes disruptive mutation, its downstream targets might still be activated by alternative routes that do not involve the disrupted element.
the observation that genes occupying more hierarchical positions in the signaling network are not subject to stronger selective forces is consistent with this emerging integrative model.
although somewhat speculative, this observation would be compatible with the organization of the signal transduction network being less hierarchical than that of biosynthetic and transcriptional regulatory networks.
in the future, a deeper understanding of the differences between the nature of biosynthetic, transcriptional regulatory, and signaling pathways may shed more light on their different patterns of evolution.
CONCLUSIONS
surprisingly, genes occupying the higher hierarchical positions of the human signal transduction network are not subject to stronger levels of purifying selection, suggesting that they are not more important for the function of the network and the fitness of the organism than genes occupying the lower hierarchical positions.
this observation sharply contrasts with the patterns observed in metabolic and transcriptional regulatory pathways and networks, in which upstream genes are generally the most selectively constrained.
these contrasting patterns of evolution might reflect fundamental differences in the function and organization of signaling and biosynthetic and transcriptional regulatory networks.
in any case, results presented here broaden our knowledge on how natural selection distributes across molecular networks.

 proteins’ rates of evolution vary across orders of magnitude. understanding the factors underlying this enormous variability is a central problem in evolutionary biology. the neutral theory of molecular evolution predicts higher selective pressures acting on genes encoding proteins that perform the most important functions for organisms’ biological fitness; however, the relative contributions of genes to fitness remain elusive. over the last decade, the availability of both genomic and functional data has allowed us to identify some of the factors that correlate with genes’ rates of evolution, among which patterns and levels of gene expression appear to be the most important. however, this gene-centered approach has been able to explain only a relatively small fraction of the variability of rates of evolution. genes and proteins do not act in isolation, but usually operate as components of complex networks of interacting molecules. therefore, taking into account the interaction networks in which a gene participates may provide key insight into the evolutionary forces acting on it. furthermore, understanding the function and evolution of molecular networks can aid applications such as metabolic engineering and drug discovery and design. molecular networks in general can be modeled as undirected graphs, whose nodes and edges represent proteins and interactions, respectively. over the last years, the accumulation of interactomic data has allowed us to cast a first glance at the impact of genes’ position in the network on their patterns of molecular evolution. for instance, genes encoding the most connected proteins in metabolic and protein–protein interaction networks are more selectively constrained than those acting at the network periphery, and interacting genes tend to show correlated evolutionary histories, e.g. evolving at similar rates. these observations indicate that the structure of molecular networks impose constraints on the evolution of their components. in addition to the centrality of proteins in functional networks, their hierarchical position in molecular pathways has also been put forward as a determinant of their rates of evolution. remarkably, in a number of biosynthetic pathways, the strength of purifying selection has been found to correlate with the position of genes along the upstream/downstream pathway axis, with upstream genes being generally the most constrained. this distribution of the levels of selective constraint has been attributed to upstream genes being involved in the biosynthesis of a broader range of biochemical compounds than downstream genes, and hence being more pleiotropic, and/or to upstream genes exerting a higher influence over the flux of metabolites along pathways. consistent observations have been made in transcriptional regulatory networks: upstream genes in these networks are more likely to be essential, and simulation analysis has shown that genes exerting a higher degree of control over other genes, and/or those that are less regulated by other genes, are more strongly affected by selection. less well understood is how selective pressures distribute along signal transduction pathways. these pathways play a key role in modulating cell function in response to extracellular and intracellular stimuli. they are characterized by the presence of a receptor that, in response to signals, is able to trigger a cascade of events that ultimately modulate the final effectors, which are responsible for mediating the biological responses of the pathway. the signal is transmitted from the upstream to the downstream part of the pathway by means of a series of molecular interactions. these interactions have a directed nature, i.e., each interaction involves an upstream modulator and a downstream modulated element. therefore, signaling pathways and networks can be represented as directed graphs, in which nodes and arcs represent molecules and directed interactions, respectively. despite the availability of some directed interaction networks, their directed nature has often been overlooked in evolutionary analyses, being often treated as undirected networks. so far, the distribution of levels of selective constraint across the upstream/downstream axis of signaling pathways has been studied in only a few pathways, with contrasting results. in the drosophila ras pathway, upstream genes appear to be more selectively constrained than those acting at the downstream part. this is to be expected under a model in which upstream genes, owing to their higher hierarchical position, have a higher impact on the function of the pathway. indeed, a single activated upstream protein has the potential to activate several downstream target molecules; therefore, higher levels of selective constraint might be expected at these genes. however, the opposite pattern has been described in the drosophila, caenorhabditis and vertebrate insulin/tor pathway, and in the yeast hog pathway. in both the yeast hog pathway and the caenorhabditis insulin/tor pathway, the tendency appears to vanish when the effect of confounding factors that correlate with both network position and levels of selective constraint, such as expression level or codon bias, are removed, whereas in the drosophila and vertebrate insulin/tor pathway the polarity in the levels of selective constraint is independent of these and other factors. therefore, the effect of upstream/downstream position of genes in signal transduction pathways on their rates of evolution, and the relative contribution of potential confounding factors, remain elusive. here, i combine comparative genomics and manually-curated directed protein interaction data to determine whether purifying selection differentially acts on the upstream and downstream genes of the human signal transduction network. for that purpose, new methods that explicitly take into account the directed character of the network are introduced. surprisingly, no evidence was found supporting that upstream genes are subject to increased levels of selective constraint. on the contrary, purifying selection was found to act more severely on genes acting at the downstream part of the network, which seems to be due to downstream genes being more highly and broadly expressed, having certain functions and, in particular, encoding proteins that are more highly connected in the protein–protein interaction network. when the effect of these confounding factors was discounted, upstream and downstream genes were found to evolve at similar rates. these results clearly indicate that the hierarchical position of genes in the signal transduction network has, in general, no direct effect on their rates of evolution, suggesting that genes occupying more hierarchical positions in the network are not more relevant for the function of the network and the fitness of the organism. the trend found in the overall signaling network is exemplified by analysis of the distribution of purifying selection along the mammalian ras signaling pathway, showing that, contrary to previously suggested in drosophila, upstream and downstream genes evolve at similar rates. these observations imply that natural selection differently distributes across signal transduction networks and across biosynthetic and transcriptional regulatory networks. this might reflect fundamental differences in the function and organization of the signal transduction network and metabolic and transcriptional regulatory networks.  unexpectedly, genes acting at the upstream part of the human signal transduction network exhibit higher rates of evolution the human signal transduction network gleaned by cui et al. was used. in this network, directed edges between pairs of proteins represent regulatory interactions — e.g., the directed edge “a → b” denotes that protein “a” activates or inhibits protein “b”. the dataset was assembled and manually curated by cui et al. by merging the interactions compiled in different databases, followed by validation of each interaction by accurate revision of the literature. undirected edges, which represent generic physical interactions, were not considered in the current analysis, since they contain no information about the hierarchical position of genes in the network. after filtering the dataset, it comprised proteins connected by non-redundant directed interactions. this dataset included ligands, receptors, adapters, kinases, phosphatases, transcription factors, structural proteins, and ion transmembrane transporters. most edges represent direct physical protein–protein interactions, with the only exception of a small fraction of “transcription factor → target gene” edges. for each gene in the network, the 1: mouse ortholog was identified and the impact of natural selection was characterized from the nonsynonymous to synonymous divergence ratio. values of ω equal to are expected for genes evolving neutrally, whereas ω <  is indicative of the action of purifying selection preserving the sequence of the encoded proteins, and ω >  in a number of codons is indicative of the action of positive selection driving the fixation of nonsynonymous substitutions. the estimated ω values range from to, indicating that purifying selection has played an important role in the evolution of these genes. i evaluated whether ω values depend on the hierarchical position that genes occupy in the network. for that purpose, three different approaches were adopted to compare the evolutionary rates of upstream and downstream network genes. first, the ω value of each gene was compared with the central value of the ω values of all its direct downstream targets. out of the genes in the network, encode proteins that have at least one direct downstream target in the network. of those, genes exhibit a higher ω value than their targets, whereas show lower ω values than their targets. this asymmetry is not compatible with a random distribution of ω values across upstream and downstream genes. therefore, surprisingly, upstream genes tend to show higher ω values than the genes that they modulate. in order to establish whether the evolution of nonsynonymous or synonymous sites is responsible for this trend, this analysis was conducted for dn and ds separately. clearly significant differences were observed for dn, but only marginally significant differences were observed for ds, indicating that the evolutionary patterns of nonsynonymous sites are the main responsible for the trend observed in ω. *, p<; **, p<, ***, p<. anumber of genes with at least one direct downstream target. only genes with available information for the parameter of interest were considered. bnumber of genes with a higher value than the central value of its downstream targets. cnumber of genes with a lower value than the central value of its downstream targets. next, evolutionary rates of genes encoding proteins occupying extreme upstream and downstream positions in the network were compared. a total of genes were classified as upstream and were classified as downstream. again, upstream genes were found to exhibit significantly higher ω values: upstream and downstream genes exhibit median ω values of and, respectively. significant differences were also observed in the dn, but not in the ds values. p-values were obtained from the mann–whitney u test. *, p < . finally, i evaluated whether ω values correlate with parameters that reflect the hierarchical position of genes in the network. both in-degree and out-degree negatively correlate with ω values, consistent with previous observations that highly connected genes are subject to stronger purifying selection. remarkably, the correlation is stronger for in-degree than for out-degree, indicating that the ω value of a gene is more influenced by the number of upstream genes, than by the number of downstream targets. similar results were obtained for dn, whereas ds does not correlate with in-degree or out-degree. beyond considering the direct downstream targets of a gene, i also considered the effect of the total number of genes acting directly or indirectly downstream of a given gene on its evolutionary rate. this parameter positively correlates with ω, again indicating that proteins potentially affecting the activity of a higher number of other proteins are encoded by genes that exhibit higher ω values. again, significant results were obtained for dn, but not for ds. *, p < ; ***, p < . taken together, these results show that evolutionary rates exhibit a polarity across the human signaling network, with downstream genes progressively exhibiting lower ω and dn values. this observation is in sharp contrast with patterns observed in metabolic and transcriptional regulatory pathways, in which upstream genes are generally the ones exhibiting lower rates of evolution. i next considered a number of potential factors that could account for the observed polarity. the higher rates of evolution of upstream genes are the result of lower levels of purifying selection rather than a higher incidence of positive selection the lower ω and dn values observed in the downstream network genes suggest a higher strength of purifying selection acting on these genes. however, ω and dn values are also affected by the action of positive selection, and hence the observed trends could alternatively be the result of a putatively higher incidence of positive selection on upstream genes. in order to discard this possibility, analyses were repeated after eliminating genes with the signature of positive selection. when applied to pairs of human-mouse orthologs, the m vs. m test identified network genes with potential signatures of positive selection. excluding these genes yielded similar results. first, after excluding these genes, genes have at least one direct downstream target in the network, out of which and show higher and lower ω values than their targets, respectively. significant results were also obtained for dn. second, the genes classified as upstream that display no signals of positive selection exhibit significantly higher ω, and dn values than the downstream genes without the signature of positive selection. finally, ω and dn values positively correlate with h, and more strongly correlate with in-degree than with out-degree. because detection of positive selection is sensitive to the number of sequences used, analyses were repeated using a total of mammalian genomes. a total of network genes exhibit a 1: orthology across all genomes, out of which exhibit potential signals of positive selection. again, excluding these genes did not alter the observations that genes acting at the downstream part of the network exhibit significantly lower ω and dn values. therefore, the distribution of ω and dn values across the network may be the result of a polarity in the strength of purifying selection, with downstream genes being subject to higher levels of selective constraint, and not the result of a higher incidence of positive selection in the upstream part of the network. subcellular location does not account for the higher rates of evolution of upstream network proteins genes acting at the extracellular compartment have been reported to exhibit high rates of evolution. consistently, among genes in the used dataset, those acting at the extracellular space exhibit higher ω and dn values than the rest of the genes in the signal transduction network. this could potentially explain the observed polarity in ω and dn values if proteins occupying highly hierarchical positions in the signaling network preferentially located to the extracellular space. however, i tested this possibility by comparing the proportion of extracellular proteins among proteins acting at the upstream and downstream parts of the network, finding no significant differences. i.e., proteins occupying extreme upstream positions in the network do not tend to act at the extracellular space, thus making it unlikely that the observed high ω values in genes encoding the upstream part of the network would be linked to subcellular location. in order to further discard this possibility, analyses were repeated eliminating the proteins targeted to the extracellular space, with similar results. first, genes are significantly less selectively constrained than their downstream direct targets. second, the non-extracellular proteins classified as upstream exhibit higher ω and dn values than the non-extracellular proteins acting at the downstream part of the network, although the differences are not statistically significant. third, the ω and dn values correlate better with in-degree than with out-degree. although evolutionary rates do not significantly correlate with h, the direction of the correlations is positive and comparable in strength to that observed in the overall dataset. these results indicate that the polarity observed in the levels of selective constraint is for the most part independent of the high evolutionary rates of proteins acting at the extracellular space. using an expanded dataset, cui, purisima and wang found that signaling genes acting at the cell membrane, cytoplasm and nucleus exhibited different rates of evolution. again, these differences raise the possibility that the observed polarity in the ω values across the network might be in part the result of the distribution of upstream and downstream proteins across the different subcellular compartments. however, i found that genes acting at the different cell compartments do not exhibit significantly different rates of evolution for the subset of genes included in the current analysis, and that the upstream and downstream parts of the network do not contain different proportions of proteins locating to the cell membrane, cytoplasm or nucleus, making it unlikely that the polarity in the levels of selective constraint described here would obey to proteins’ subcellular location. in order to further discard this possibility, analyses were conducted separately on subsets of genes acting at the different cell locations, with consistent results across all four compartments. first, genes tend to evolve faster than their direct downstream targets acting at the same cell compartment, when the analyses are restricted to genes acting at the extracellular space, the cell membrane, the cytoplasm, or the nucleus. second, ω values are higher for upstream than for downstream genes when genes acting at the extracellular space, the cell membrane, the cytoplasm or the nucleus are analyzed separately, although the differences are not statistically significant. taken together, these results indicate that the polarity in the levels of selective constraint observed across the network is independent of the subcellular location of its components. the lack of significance in some tests probably owes to the reduced statistical power resulting from partitioning the dataset. gene functions upstream and downstream of the signal transduction network genes with different functions are known to evolve at different rates. hence, in order to gain further insight into the different selective pressures acting at the upstream and downstream parts of the human signaling network, it was tested whether the genes acting at the upstream part performed different functions to the genes acting at the downstream part. for that purpose, the frequencies of the goslim goa terms of both groups was compared using the fatigo software. upstream genes exhibited a significant enrichment in the term “kinase activity”, whereas downstream genes were enriched in the category “ion transmembrane transporter activity”. hence, the higher levels of selective constraint observed in the downstream part of the signaling network could potentially result from this unequal distribution of gene functions, if genes encoding kinases evolved particularly fast, and/or those encoding ion transmembrane transporters evolved particularly slow. both gene groups were found to evolve slower than the rest of genes in the network, raising the possibility that the higher levels of selective constraint acting at the downstream part of the network would be due to the enrichment of this part of the network in ion transmembrane transporters. in order to discard this possibility, analyses were repeated excluding the network genes encoding ion transmembrane transporters, with similar results, although not statistically significant in all cases. first, genes are significantly less selectively constrained than their downstream direct targets. second, the non-ion transporter proteins classified as upstream exhibit higher ω and dn values than the non-ion transporter proteins acting at the downstream part of the network, although the differences are not statistically significant. third, the ω and dn values correlate better with in-degree, than with out-degree, and they exhibit a positive, although non-significant, correlation with h. when genes encoding kinases were also excluded from the analysis, the correlation between h and rates of evolution reached significance. these results indicate that, although the different functions of upstream and downstream genes may partially contribute to the observed polarity in the levels of selective constraint, this factor does not completely account for the trend. genes acting at the downstream part of the network exhibit higher expression levels and breadths, and encode more highly connected proteins levels of selective constraint acting on a gene are affected by a number of factors, including gene expression level and breadth, codon bias, connectivity of the encoded proteins, and the number of paralogs. measures of these parameters were obtained from different sources. because the dataset generated by cui et al. focuses on signaling molecules, proteins’ total connectivities were measured as the total number of physical interactors in the whole human protein–protein interaction network. all these factors significantly correlate with ω and dn in the dataset used. therefore, a putative polarity in the distribution of these factors across the network could potentially account for the observed polarity in the levels of selective constraint. remarkably, expression breadth and connectivity seem to be the best correlates of the rates of evolution, with both variables exhibiting a similar degree of correlation with ω and dn. this result contrasts with previous observations that connectivity is a relatively weak predictor of rates of evolution, suggesting that this correlation is particularly strong for genes involved in signal transduction. furthermore, the correlation between connectivity and rates of evolution remains significant when expression level and breadth, codon bias, and number of paralogs are simultaneously controlled for. although genome-level analyses have revealed that the length of the encoded proteins also correlates with ω and dn, these correlations are not significant in the gene set used, and hence this parameter was not considered. **, p<, ***, p<. i evaluated whether factors correlating with ω and dn differ between upstream and downstream genes. the paired test shows that genes tend to exhibit a lower expression breadth, and to encode proteins with a lower connectivity than their direct downstream targets. the test is only marginally significant for expression level. as expected, both in-degree and out-degree positively correlate with connectivity. the correlation is stronger for in-degree than for out-degree, consistent with a higher connectivity in the downstream part of the network. genes occupying extreme downstream positions in the pathway exhibit higher levels and breadths of expression, and encode more highly connected proteins, than upstream genes, although none of these differences is statistically significant. the lack of significance in the mann–whitney tests might result from a reduced statistical power resulting from the fact that only a fraction of the genes were used in this analysis, whereas all genes were used in the other two analyses. although none of the studied factors significantly depends on the upstream/downstream position of genes in all three analyses, these results provide evidence that expression breadth and connectivity, and perhaps expression level, exhibit a polarity across the human signaling network, with downstream genes being more highly and broadly expressed, and encoding more highly connected proteins than upstream genes. highly and broadly expressed genes, and those encoding highly connected proteins, tend to be highly constrained. combined, these trends raise the possibility that the polarity observed in the levels of selective constraint could be the result of the distribution of these factors across the network. the lower evolutionary rates of downstream network genes are a byproduct of the distribution of protein connectivity and gene expression level and breadth across the network i considered whether the observed polarity in the levels of selective constraint across the human signaling network was the result of the distribution of factors correlating with ω and dn. for that purpose, this polarity was evaluated after discounting the putative effects of these factors. linear regression was used to model the dependence of ω and dn from each factor separately, and the residuals of the model for each gene were used in the paired and the mann–whitney tests. the paired tests show that genes exhibit significantly higher ω and dn values than their downstream targets even if the effect of expression level, expression breadth, codon bias, or the number of paralogs is discounted, indicating that the difference in the intensity of purifying selection acting on upstream and downstream genes is independent of these factors. however, when the effect of connectivity is discounted, the test yields no significant results. genes occupying extreme upstream positions in the network are significantly less constrained than those acting at the downstream part, even when the effect of expression breadth, codon bias, or the number of paralogs is removed, but the difference is not significant when expression level or connectivity are included in the analysis. given that the relationship between evolutionary rates and their correlates is not necessarily linear in all cases, a complementary nonparametric analysis was performed to assess the differences in ω and dn values between genes occupying extreme upstream and downstream positions. partial correlation analysis was used to evaluate the association between upstream/downstream position and rates of evolution while controlling for each of the controlling variables. similar results were obtained: upstream and downstream genes exhibit significantly higher ω and dn values when controlling for codon bias or the number of paralogs, but not when controlling for expression level, expression breadth, or connectivity. finally, partial correlation analysis was used to evaluate the association between h and ω and dn while controlling separately for each factor. both correlations are significant when codon bias or the number of paralogs are used as controlling variables, but not when expression level, expression breadth, or connectivity are included in the analysis. therefore, the differences in the ω and dn values of upstream and downstream network genes disappear when the effect of expression level, expression breadth, or connectivity, are removed. this indicates that these factors account for the higher levels of selective constraint acting on the downstream genes of the human signal transduction network. in particular, when connectivities are controlled for, ω and dn values are not significantly different for upstream and downstream genes, regardless of the method used, indicating that this factor may be the main responsible for the polarity in the levels of selective constraint. when the effects of expression level or breadth are accounted for, the differences in ω and dn between upstream and downstream genes vanish for only one or two of the three methods used, indicating that these factors may contribute to the selective constraint polarity to a lesser extent. a case study: natural selection upstream and downstream of the mammalian ras signaling pathway so far, the drosophila ras pathway is the only signal transduction pathway for which higher levels of selective constraint have been described in genes acting at the upstream part. the ras pathway is activated by diverse extracellular stimuli, including growth factors, cytokines and hormones. activation of receptor tyrosine kinases triggers a cascade of signaling events that ultimately promote cell proliferation, differentiation, migration and survival. the architecture of this pathway has been characterized in detail in a number of organisms, including drosophila and mammals, revealing a highly conserved structure and function across metazoans. riley et al. examined the patterns of molecular evolution of six of the genes involved in the drosophila ras pathway. for that purpose, they surveyed the levels and patterns of polymorphism in drosophila melanogaster, and divergence with respect to the sister species d. simulans. they found that three proteins acting at the upstream part of the pathway were identical in both species, and exhibited no or just a few rare replacement polymorphisms within d. melanogaster. on the contrary, the three downstream proteins included in their analysis displayed a number of fixed differences between the two species, as well as some replacement polymorphisms in d. melanogaster. consistent results were obtained from the ω values, pointing out to higher levels of purifying selection acting at the upstream part of the pathway. a proposed explanation for this observation was that, owing to their highly hierarchical position, proteins acting at the upstream part of the pathway would act as key control points, greatly influencing the activity of downstream proteins, whereas downstream proteins would act as modulators of the output signal. the current wealth of genomic data provides an excellent opportunity to study the patterns of molecular evolution of all the core components of the ras pathway in other taxonomic groups. hence, i considered whether levels of selective constraint acting on the mammalian ras pathway genes exhibit a similar distribution across the pathway to that observed in drosophila. most of the afore-mentioned drosophila pathway genes have multiple orthologs in mammalian genomes known to participate in the ras pathway. for each such mammalian gene, the ω value was estimated from comparison of the human-mouse pair of orthologs. genes encoding the mammalian proteins sos sos erk erk and rsk1– which are also key components of the pathway, were also included in the current analysis. contrary to previously suggested in drosophila, genes acting at the upstream part of the mammalian ras pathway were not found to exhibit lower rates of evolution than those acting at the downstream part. if anything, the opposite seemed to be the case, although the differences between both groups were not significant. consistently, ω values were found to exhibit a negative, although non-significant, correlation with the position of genes along the upstream/downstream axis of the pathway. for this analysis, pathway position was defined as the number of steps required to transduce the signal from the most upstream element, to each of the pathway components. 1: table s for a full list of the core pathway genes and their rates of evolution. neither expression levels, nor expression breadths or connectivities significantly correlate with pathway position or with ω values. it should be noted, however, that microarray-based gene expression data, such as the one used here, is subject to high levels of measurement noise. more accurate measures of expression are available for some of the proteins involved in the ras pathway, and indeed it has been noticed that downstream proteins are present in higher concentrations in the cell, thus leaving open the possibility that the slightly higher levels of selective constraint observed in the downstream part of the pathway may result from these genes being expressed at higher levels. in any case, these results indicate that, contrary to previously suggested in drosophila, genes acting at the upstream part of the mammalian ras pathway are not subject to higher levels of selective constraint than those acting at the downstream part. if anything, the opposite might be the case, although the trend is not statistically significant. this observation is in agreement with the general trend described here for the whole human signal transduction network. discussion in summary, results presented here demonstrate that genes acting at the upstream part of the human signal transduction network clearly do not evolve under stronger selective pressures. this finding seems striking, as these proteins, owing to their higher hierarchical position, may be capable of regulating the activity of a higher number of downstream proteins. therefore, mutations in the genes encoding upstream network proteins would be expected to have more pleiotropic effects, and hence these genes would be expected to evolve under stronger purifying selection. in contrast with these expectations, a number of analyses that capture different aspects of the hierarchical position of genes in the signaling network found no evidence for upstream genes being subject to higher levels of selective constraint. on the contrary, upstream network genes are less selectively constrained than downstream genes, which appears to be explained by downstream genes performing certain functions, being more highly and broadly expressed and, in particular, encoding more highly connected proteins. protein connectivity, gene expression breadth and, to a lesser extent, expression level, exhibit a polarity across the human signaling network, with downstream genes progressively encoding proteins that are more connected, and being expressed in a higher number of tissues and at higher levels. the higher levels of expression of downstream genes might facilitate the amplification of signals through their progression throughout signaling cascades, which requires higher protein abundances for the downstream elements of the pathway. additionally, the higher connectivity of downstream genes is consistent with previous observations that hubs tend to be tightly regulated, exhibiting higher numbers of phosphorylation sites and higher mrna decay rates. the distribution of these factors across the network, combined with the fact that highly and broadly expressed genes, and those encoding highly connected proteins in the protein–protein interaction network, tend to evolve under stronger purifying selection, accounts for the higher levels of selective constraint observed in genes acting at the downstream part of the network. in particular, when the effect of connectivity is factored out, the differences in ω and dn values between upstream and downstream genes are no longer significant, regardless of the method used, pointing to this factor as the main responsible for the higher levels of selective constraint acting on the downstream network genes. when the effects of expression level and breadth are removed, the differences in the levels of selective constraint observed between upstream and downstream genes also vanish for some of the tests, suggesting that these factors also contribute to the tendency. likewise, protein functions may also play a role, as the trend also vanishes for some tests when ion transmembrane transporters are removed from the analyses. in any case, taken together these observations indicate that the upstream/downstream position of proteins in the human signal transduction network does not have a direct effect on their rates of evolution. the global trends observed in the entire human signal transduction network mirror those previously observed in some particular pathways. in both the caenorhabditis insulin/tor signaling pathway and the yeast hog signaling pathway, downstream genes evolve under stronger levels of selective constraint, probably as a result of their higher levels of expression and codon bias, respectively. it should be noted, however, that the general trends observed across the whole network do not necessarily imply that each particular signaling pathway must exhibit the same trend. for instance, in the drosophila and vertebrate insulin/tor pathways, downstream genes evolve under stronger purifying selection, even when the putative effect of several factors that correlate with levels of selective constraint is discounted. another example that does not seem to conform to the general trend described here is the drosophila ras pathway, whose upstream genes appear to be the most constrained, despite these genes being the ones expressed at lower levels. therefore, despite the overall trends described here, the particular structure and function of each individual pathway may result in different distributions of the impact of evolutionary forces. for instance, pathways whose upstream genes, for some reason, present higher connectivities, and/or levels or breadths of expression, would probably be more selective constrained in the upstream part. it should also be noted that factors other than the ones considered here may also affect the distribution of levels of selective constraint across functional pathways and networks. in particular, enzymes whose kinetic characteristics exert a high degree of influence over the entire network are expected to be more selectively constrained than those exerting little influence on the network dynamics. the degree to which a particular enzyme can influence the behavior of the entire network is, however, difficult to determine from currently available data. so far, the drosophila ras pathway is, to the best of the author’s knowledge, the only signal transduction pathway whose upstream genes appear to be the ones evolving under higher levels of purifying selection. however, due to the limited availability of genomic data at that time, this observation was based on the analysis of a small set of genes and the comparison of two closely related species. results reported here for all core components of the mammalian ras pathway, and based on human-mouse comparisons, do not point in the same direction: upstream genes do not exhibit lower rates of evolution than those acting at the downstream part of the pathway. if anything, the opposite might be the case. in any case, this observation lends further support to the finding that genes occupying highly hierarchical positions in the human signal transduction network are not subject to higher levels of selective constraint than those acting at the downstream part. the observation that genes acting at the upstream part of the human signal transduction network are not subject to higher levels of selective constraint sharply contrasts with the trends generally observed in biosynthetic pathways. these pathways consist of a series of enzymes that share metabolites, with the product of a given enzyme being usually the substrate of the following enzyme in the pathway. in such pathways, metabolites can often take different alternative routes, which end up in the biosynthesis of different end products. in a growing number of biosynthetic pathways, including the plant anthocyanin, isoprene, terpenoid, and carotenoid pathways, the bombyx melanin pathway, and the mammalian dopamine pathway, levels of selective constraint correlate with the position of genes along the upstream/downstream pathway axis, with upstream genes being always the most constrained. two plausible scenarios might explain such a distribution of selective constraints in biosynthetic pathways. first, as a result of their branching topology, these pathways often exhibit a hierarchical structure, with upstream genes being necessary for the biosynthesis of a higher number of end products than downstream genes; therefore, mutations in genes acting at the upstream part may have higher pleiotropic effects. second, upstream enzymes in biosynthetic pathways probably often exert a higher control over the flux of metabolites throughout the pathway; thus, the overall pathway function may more strongly depend on the function of these enzymes. both models predict a higher relevance for pathway function of upstream genes, and hence a higher strength of natural selection acting on these genes. similar trends have been observed in transcriptional regulatory networks. first, genes occupying higher hierarchical positions in these networks are likely to be essential, reflecting their key importance for organism’s fitness. second, simulation analyses of transcriptional regulatory networks have shown that genes exerting a high degree of control over other genes, or those that are less regulated by other genes, are expected to be more strongly affected by selection. therefore, results presented here imply that the distribution of evolutionary forces across the upstream and downstream parts of the signal transduction network remarkably differs from that observed in biosynthetic and transcriptional regulatory pathways. in biosynthetic and transcriptional regulatory pathways, upstream genes evolve under stronger levels of selective constraint, whereas this is not the case for the human signal transduction network. this suggests that genes occupying higher hierarchical positions are not more relevant for the function of the signaling network, or to the organism’s fitness, than downstream genes. these different patterns of evolution may reflect fundamental differences in the function and organization of the different kinds of pathways. as signaling networks become better understood, it becomes more apparent that signal transduction does not take place in the context of independent linear pathways in which each protein is exclusively controlled by the preceding protein in the cascade. this classical model is being progressively substituted by a much more integrative model of signal transduction, in which a considerable amount of cross-talk takes place between pathways. as a result, pathways would cooperate to integrate signals, and to make decisions about the responses in a decentralized manner, thus conferring robustness to the system. this kind of organization may allow the network to accommodate modifications in its upstream genes. for instance, even if an upstream gene undergoes disruptive mutation, its downstream targets might still be activated by alternative routes that do not involve the disrupted element. the observation that genes occupying more hierarchical positions in the signaling network are not subject to stronger selective forces is consistent with this emerging integrative model. although somewhat speculative, this observation would be compatible with the organization of the signal transduction network being less hierarchical than that of biosynthetic and transcriptional regulatory networks. in the future, a deeper understanding of the differences between the nature of biosynthetic, transcriptional regulatory, and signaling pathways may shed more light on their different patterns of evolution.  surprisingly, genes occupying the higher hierarchical positions of the human signal transduction network are not subject to stronger levels of purifying selection, suggesting that they are not more important for the function of the network and the fitness of the organism than genes occupying the lower hierarchical positions. this observation sharply contrasts with the patterns observed in metabolic and transcriptional regulatory pathways and networks, in which upstream genes are generally the most selectively constrained. these contrasting patterns of evolution might reflect fundamental differences in the function and organization of signaling and biosynthetic and transcriptional regulatory networks. in any case, results presented here broaden our knowledge on how natural selection distributes across molecular networks. 
these results indicate that the upstream/downstream position of proteins in the signal transduction network has, in general, no direct effect on their rates of evolution, suggesting that upstream and downstream genes are similarly important for the function of the network. in general, proteins whose function has a great impact on organisms’ fitness are expected to evolve under stronger selective pressures. less well understood is, however, how selective pressures distribute along signal transduction pathways. an open question in evolutionary biology is what factors are responsible for this variability. here, i combine comparative genomics and directed protein interaction data to study the distribution of evolutionary forces across the human signal transduction network. this implies that natural selection differently distributes across signal transduction networks and across biosynthetic and transcriptional regulatory networks, which might reflect fundamental differences in their function and organization. the trends found in the overall signaling network are exemplified by analysis of the distribution of purifying selection along the mammalian ras signaling pathway, showing that upstream and downstream genes evolve at similar rates. in biosynthetic pathways, upstream genes usually evolve under higher levels of selective constraint than those acting at the downstream part, as a result of their higher hierarchical position. when the effect of these confounding factors is discounted, upstream and downstream genes evolve at similar rates. on the contrary, purifying selection was found to act more strongly on genes acting at the downstream part of the network, which seems to be due to downstream genes being more highly and broadly expressed, performing certain functions and, in particular, encoding proteins that are more highly connected in the protein–protein interaction network. proteins evolve at disparate rates, as a result of the action of different types and strengths of evolutionary forces. surprisingly, no evidence was found for higher levels of selective constraint at the upstream network genes. network evolutionevolutionary rateselective constraintpurifying selectionsignal transductionmolecular pathwaysupstream and downstream position similar observations have been made in transcriptional regulatory networks, whose upstream elements appear to be more essential and subject to selection.
 proteins evolve at disparate rates, as a result of the action of different types and strengths of evolutionary forces. an open question in evolutionary biology is what factors are responsible for this variability. in general, proteins whose function has a great impact on organisms’ fitness are expected to evolve under stronger selective pressures. in biosynthetic pathways, upstream genes usually evolve under higher levels of selective constraint than those acting at the downstream part, as a result of their higher hierarchical position. similar observations have been made in transcriptional regulatory networks, whose upstream elements appear to be more essential and subject to selection. less well understood is, however, how selective pressures distribute along signal transduction pathways.  here, i combine comparative genomics and directed protein interaction data to study the distribution of evolutionary forces across the human signal transduction network. surprisingly, no evidence was found for higher levels of selective constraint at the upstream network genes. on the contrary, purifying selection was found to act more strongly on genes acting at the downstream part of the network, which seems to be due to downstream genes being more highly and broadly expressed, performing certain functions and, in particular, encoding proteins that are more highly connected in the protein–protein interaction network. when the effect of these confounding factors is discounted, upstream and downstream genes evolve at similar rates. the trends found in the overall signaling network are exemplified by analysis of the distribution of purifying selection along the mammalian ras signaling pathway, showing that upstream and downstream genes evolve at similar rates.  these results indicate that the upstream/downstream position of proteins in the signal transduction network has, in general, no direct effect on their rates of evolution, suggesting that upstream and downstream genes are similarly important for the function of the network. this implies that natural selection differently distributes across signal transduction networks and across biosynthetic and transcriptional regulatory networks, which might reflect fundamental differences in their function and organization. network evolutionevolutionary rateselective constraintpurifying selectionsignal transductionmolecular pathwaysupstream and downstream position 