The function of communities in protein interaction networks at multiple scales

Protein-Protein Interaction Datasets

Here we use the BioGrid (http://www.thebiogrid.org, downloaded January 2010, [37]), IntAct (http://www.ebi.ac.uk/intact, downloaded January 2010, [38]) and Mint databases (http://mint.bio.uniroma2.it/mint, downloaded January 2010, [39]) to assemble our protein interaction networks. We use only interactions between proteins that have an SGD identification (Saccharomyces Genome Database, http://www.yeastgenome.org, [40]).

We divide interactions on the basis of their type (A or P) and hence assemble the two networks. The IntAct database [38] gives interaction types from the Molecular Interaction ontology [33] directly. It contains 23632 interactions of type A and 26611 of type P. The Mint database [39] uses the Molecular Interaction interaction detection type ontology, the broad categories of which are biophysical, biochemical, and protein complementation assay. The biochemical techniques give evidence of association (type A interactions), and the biophysical and protein complementation assays give evidence of physical interactions (type P). Using this division, there are 13347 A type interactions and 10407 P type interactions. The BioGrid database [37] uses its own evidence types. Those giving evidence of P type interactions are reconstituted complex, PCA, Co-crystal structure and yeast-two-hybrid. Those giving evidence of type A interactions are affinity capture, biochemical activity, co-fractionation, co-purification and Far Western. (Details of these experimental types can be found on the BioGrid website, http://www.thebiogrid.org). There are 35716 A type interactions and 13142 P type interactions overall. Of the potential 6607 proteins in the yeast proteome http://www.yeastgenome.org, there are 5002 proteins connected by A type interactions, and 5692 connected by P type interactions. Here we only study the largest connected component of these networks, leaving 4980 proteins in the A network and 5669 in the P network. Some summary statistics for the two amalgamated networks are shown in Table 1. The A network is denser, and has higher clustering. There are 5947 interactions in common between the A and the P networks.

Potts community detection

We apply the Potts method [23]. It partitions the proteins into communities at many different values of a resolution parameter, thus finding communities at different scales within the network. The method seeks a partition of nodes into communities that minimises a quality function ('energy'):

(1)

where s _i is the community of node i, δ is the Kronecker delta, λ is the resolution parameter, and the interaction matrix J _ij (λ) gives an indication of how much more connected two nodes are than one would expect at random (i.e., in comparison to some null hypothesis). The energy H is thus given by a sum of elements of J for which the two nodes are in the same community. Optimising H is known to be an NP-hard problem [41, 42], so one must use a computational heuristic. Here we use the greedy algorithm discussed in [43] and freely available http://www.lambiotte.be/codes.html, which performs well against various benchmark tests [44]. As pointed out by Good et al [45], one must be cautious in interpreting results obtained from detecting communities by optimizing modularity or similar quality functions, as there is a degeneracy of partitions with almost optimal H. As a consequence different optimisation techniques can find very different optima.

The interaction matrix J has elements

(2)

where the matrix B with elements B _ij is the adjacency matrix. In this case B _ij = 1 if proteins i and j interact, and B _ij = 0 otherwise. The matrix R with elements R _ij defines a null model, against which we are comparing the network of interest. Here we choose the standard Newman-Girvan null model [46], which has the property that it preserves the expected node degree sequence. That is,

(3)

where k _i = ∑ _j B _ij is the degree of node i, and W =∑ _ij B _ij /2 is the number of edges in the network. When λ = 1, H is the standard Newman-Girvan modularity quality function, upon which many community detection algorithms are based [11, 46]. We hence refer to this value of the resolution parameter as the standard resolution. Values of λ > 1 probe the network at resolutions above the resolution limit.

We investigate partitions of the network in the range 0.1 is ≤ λ ≤ 1000, and sample at intervals of 0.01 on a logarithmic scale (we hence report results for -1 ≤ log(λ) ≤ 3). At λ = 0, all nodes in our set will be assigned to the same community. As we increase λ, communities split and become smaller. If we allow λ to increase until all of the entries in J _ij are negative, then each node will be assigned to its own community.

Convention for identifying communities at different partitions

To relate the partition at one value of the resolution parameter λ to that at another (which we use here for visualisation), we require a convention for labelling communities. Here we use a method based on the overlap of shared nodes [47]. A convention based on links rather than nodes gives nearly identical results. Let the communities in the first partition (which here is that at the highest resolution) be labeled K ₁ , ..,K _s , and those in the next partition be labelled L ₁ , ...,L _t . Then for each pair of communities, {K _i , L _j }, we have

(4)

where |B| denotes the cardinality (number of elements) of the set B. Starting with the largest value of W _ij , we relabel community i as community j. Relabelling proceeds with the next largest W _ij , as long as community i is not yet relabelled, until all communities have been relabelled. If s > t, we introduce a new label.

Pairwise measures of functional similarity

It is impossible to uniquely quantify similarity in biological function. Here we rely primarily on the GO http://www.geneontology.org, which provides the most comprehensive available database of functional annotations. We use the Biological Process sub-ontology annotations to yeast, which are maintained by the SGD consortium [40]. Terms are related to each other through a directed acyclic graph (DAG). Proteins are annotated with the most specific terms that are known about them. It is then possible to add to this set their parent terms by following the structure of the DAG, up to the root node. Well-characterised proteins are those annotated with terms far from the root node. Of the 6346 yeast proteins in the GO annotation set, 5347 have biological process annotations (excluding the root node). We carried out the same tests using the Molecular Function and Cellular Component sub-ontologies, which gave similar results.

We also use MIPS terms (http://www.helmholtz-muenchen.de/en/ibis, [21]), which are a useful double check on our results from GO, and have the added advantage that the terms are all found at the same level within the hierarchy of terms. Here we only use the top level of the MIPS hierarchy.

Following [48], we quantify the functional similarity between two proteins i and j by finding the set of GO terms annotated to both proteins and counting the total number of proteins, n _ij , that share that set of terms. We then define a similarity measure between proteins i and j as

(5)

where N is the total number of proteins. If both proteins are annotated with a set of terms that few proteins share, then they will be judged as functionally similar under this measure. Unlike many other measures, G _ij does not penalise proteins for lack of annotation when judging their similarity. This is desirable, as we know that the GO annotations (even for the well-characterised S. cerevisiae) are far from complete. The quantity M _ij is similarly defined through Equation 5 for the MIPS annotations.

The benefit of using a pairwise similarity measure that takes into account the full set of functional information available, rather than examining enrichment of function on a term by term basis, is that the measure has the potential to capture more general functional similarities between a pair of proteins.

We also define a similarity between two proteins from a single high-throughput experiment via the growth rates of knock-out strains under a range of different conditions. Using the data in [32], we define C _ij , the correlation in growth rates of the strain with gene i knocked out to the strain with gene j knocked out under 418 different conditions:

(6)

where the elements of the vector L _i are

(7)

the parameter is the mean growth rate of strain i under different control conditions, and is the growth rate under one of the 418 treatment conditions. We use the results from the homozygous strains. Because many gene deletions are lethal, there is only data available for 3625 proteins, of which 3184 are in the A network and 3422 are in the P network.

Assessment of a community's functional homogeneity

As mentioned previously, a fair test of the functional homogeneity of a community must take into account the fact that a pair of proteins that interact will be more similar than a randomly chosen pair. Standard enrichment tests do not take this into account, as they compare enrichment in a group of proteins, in this case a community, to what one would expect to attain from a randomly chosen set of proteins [49]. A community necessarily contains many more interacting pairs than a randomly chosen set. We thus compare the pairwise functional similarities of all interacting pairs of proteins in a community to the same measure for all interacting pairs in the network, thereby controlling for the number of interacting pairs.

To capture the pairwise similarity between two proteins that interact {ij}, we use z-scores:

(8)

Where S stands for one of our three similarity measures (based on GO, G, MIPS, M, or correlated growth rates, C), μ is the mean and α the standard deviation of all of the elements of S for which proteins i and j interact in the network of interest (A or P).

A desirable quality for our test of functional homogeneity is the ability to compare communities found at different resolutions in an even handed manner. It is inherent in the nature of a statistical test that the significance of the test statistic under consideration (for example, the difference between the sample mean and the population mean) depends on the sample size: if one has a larger sample size, one can judge smaller differences to be 'significant'. To determine the aggregate z-score, z_agg, for the mean of a set of individual z-scores, z_ind, one calculates , where N is the number of z_inds and μ(z_ind) is their mean [50]. So, given a μ(z_ind), a larger and hence more significant z_agg is achieved for a larger sample size (i.e. larger N). In order to separate out the effects of the number of interactors in the community from functional homogeneity, we thus choose to base assessment of functional homogeneity on the μ(z_ind), in our case μ(z_{ij}) (z_{ij}is defined in Equation 8). We judge as 'significant' all those communities that have μ(z_{ij}) above 0.3, and call such communities "functionally homogeneous". We stress that this is not strictly an assessment of statistical significance, as we are choosing to ignore sample size. The value of 0.3 would be judged to be significant at the 0.05 significance level for any community with 30 or more interacting pairs.

Classification of protein types

We focus on a small but broad set of protein types, which are the GO biological process terms within the yeast GO slim [51] that are annotated to at least 200 yeast proteins. They are (numbers of proteins in brackets): 1. DNA metabolic process (357); 2. protein modification process (465); 3. transport (859); 4. response to stress (458); 5. membrane organization (208); 6. RNA metabolic process (715); 7. vesicle-mediated transport (280); 8. response to chemical stimulus (298); 9. cellular lipid metabolic process (204); 10. cellular carbohydrate metabolic process (220) and 11. chromosome organization (338).

Topological properties that correlate well with functional homogeneity

We investigate 26 topological properties of the identified communities and assess whether any of these can be used to identify functionally homogeneous communities. Examples include mean clustering coefficient, betweenness measures, and network diameter. Any topological properties that correlate well with functional homogeneity can then be used to predict functionally homogeneous communities. We use each topological property as a classifier by predicting communities as functionally homogeneous when the value of that property is above a threshold, which we vary to construct a Receiver Operating Characteristic (ROC) curve. An ROC curve plots the number of communities correctly predicted as functionally homogeneous versus the number falsely predicted [52]. We calculate the area under the ROC curve (AUC) for each metric at each value of λ, and report the mean of this quantity over resolutions between 0 ≤ log(λ) ≤ 3 (we exclude -1 ≤ log(λ) < 0, as the results are very noisy due to the small number of communities present). An AUC of 0.5 would be expected from a random classifier. AUCs of greater than 0.5 imply that higher values of the metric are predictive of functional homogeneity. AUCs of less than 0.5 imply predictive power if below a threshold of that particular property was used (i.e. that the property and functional homogeneity are negatively correlated).

Pairwise properties of proteins

Communities

Figure 1 shows the communities that we find in the A and P yeast networks as the resolution parameter λ is varied. As λ increases, more and smaller communities are found (see Table 3). At λ = 1 (i.e. log(λ) = 0), which corresponds to standard Newman-Girvan modularity [46], most communities contain a few hundred proteins. By log(λ) = 3 however, almost all proteins are in communities of size three or smaller. As shown in Figure 1, some sets of nodes are classified in the same community through large changes in the resolution parameter and hence represent particularly inter-connected parts of the network. Figure 1 can be contrasted with Figures S1 in Additional File 1, which are similar calculations on a random network and a network designed to possess strong communities. In the former, not much structure is present, in the latter, there are very distinct blocks.

log(λ)	mean size of communities
	A	P
-0.5	681	2834
0	293	405
0.5	73	79
1	22	26
1.5	11	10
2	6	6
2.5	5	5
3	4	4

The black lines in Figure 2 illustrate for a) the A network and b) the P network, i) the number of communities of size four or more as the resolution changes, and ii) how many proteins are in those communities.

The two networks, A and P, contain very different types of interactions, and they can therefore be used to identify different aspects of the cell's functional organisation. The A network is also much denser than the P network. A interactions would therefore dominate the clustering into communities, thereby making it very hard to pick out any structures given by P type interactions (as occurs in [53]). When considering a particular protein or set of proteins, comparisons between communities found in the A and P networks can be made, see the Examples section. Global comparisons between the partitions of the A and P networks at a particular resolution are not necessarily meaningful as, for example, the size of communities depends both on the size and other properties of the network.

Data files containing the A and P networks and the community membership of proteins at multiple resolutions are available at http://www.stats.ox.ac.uk/research/proteins/resources.

Functional homogeneity of communities

We now assess how many communities are judged functionally homogeneous, looking in particular at how our results vary with resolution parameter.

Figures 2i) illustrate the number of communities judged to be functionally homogeneous, and Figures 2ii) show the number of proteins in communities judged to be functionally homogeneous, for a) the A network and b) the P network. We find that the large communities present at small values of the resolution parameter λ are not judged to be functionally homogeneous. As λ is increased, larger numbers of proteins occur in functionally homogeneous communities, peaking in the range 1.5 < log(λ) < 2. At log(λ) = 1.5, the mean community size is 73 proteins, and the majority of proteins, 3071 of 4980, are in functionally homogeneous communities as judged by our GO similarity measure. The shapes of the curves of both Figure 2a) and 2b) for all three similarity measures are very similar. Indeed, we find that the overlap between the communities judged to be functionally homogeneous between any two of the three measures is high; for example, it is 70% between the GO and correlated growth rates measure over almost the entire range of the resolution parameter in both A and P networks (see Figure S2 in Additional File 1 for the complete data). Given that the correlated growth similarity measure represents a very different data type to the GO and MIPS annotations, this agreement gives us confidence in the similarity measure we use for GO and MIPS. As we use only the top level of the MIPS functional annotations, we capture less information than the GO measure, so it is unsurprising that fewer communities are found to be functionally homogeneous under this measure.

The P network shows a similar pattern to the A network. One difference is that communities start to be judged as functionally similar at a slightly lower resolution. This is most likely due to the different topological properties of the P network. That there are comparably many functionally homogeneous communities in the P network as the A network is of interest, as communities found in P networks are found to be poor choices for predicting function on the basis of enrichment of terms [31].

For almost all proteins, there is some value of the resolution parameter that assigns them to a functionally homogeneous community. In fact 4652 out of 4980 A proteins and 5647 and of 5669 P proteins are in such communities at some value of the resolution parameter. For a given protein, it may not be that it interacts most closely with proteins involved in the same process. Indeed it is often necessary to look at a larger scale, placing the community in a bigger community in order to identify the biological processes it participates in. Whether or not this is the case, and which network scale (resolution) is most indicative of the processes a protein is involved in, will depend on the particular protein one is interested in. This demonstrates the biological motivation for investigating community structure at multiple resolutions, and suggests the desirability of a method to easily identify those communities most likely to be functionally homogeneous.

We might expect proteins involved in particular processes to show different propensities to lie in functionally homogeneous communities. We focus on a set of general protein types (as defined and listed in the Methods), and investigate what fraction of each type of protein lie in communities judged functionally homogeneous under the GO measure through changing resolution parameter. Figure 3 illustrates for the A network these percentages for four particular processes. (Figure S3 in Additional File 1 shows the same figure for all 11 terms for the A network and separately for the P network). Proteins of some types are far more likely to be found in functionally homogeneous communities than others. For example, for both the A and P networks, proteins involved in chromosome organisation are far more likely to be found in functionally homogeneous communities than proteins involved in lipid metabolism. In addition, there are some indications that the resolutions of most interest can depend on the type of protein under investigation. As can be seen in Figure 3, proteins involved in RNA metabolic processes are more likely to be found in functionally homogeneous communities at log(λ) = 0.8, where the mean size of communities is 30. In contrast, proteins involved in vesicle-mediated transport are found in greater numbers in functionally homogeneous communities at log(λ) = 1.7, where the mean size of communities is 10.

Examples of communities found at multiple resolutions

Consider the community at log(λ) = 0 that is marked as the blue block in Figure 1 for the A network (over node labels approximately 0 to 500). This contains 528 proteins and consists largely of proteins with some relationship to the ribosome (based on short protein descriptions found on the SGD website). Figure 4a) shows this community, where we have coloured nodes according to the community partition at the later partition log(λ) = 0.5. The colours - red, yellow, and blue - are the same as in Figure 1, where most of the community present at log(λ) = 0 has split into three communities at log(λ) = 0.5. The blue community consists of 107 proteins, which are largely precursors to and processors of the large ribosomal unit. The red community consists of 95 proteins, which have a similar function but for the small ribosomal subunit. The yellow community has 190 proteins, 93 of which are constituents of the ribosome and the remainder of which are either of unknown function or associate to the ribosome. We give short descriptions of the proteins in these communities in Additional File 2.

An illustration of the biological relevance of community structure at three partitions is given in Figures 4b) and 4c). We show a community of 90 proteins at log(λ) = 0.5, and display its partition into communities at b) log(λ) = 0.75 and c) log(λ) = 1.6. Almost all of the proteins in the community at log(λ) = 0.5 play some role in transcription initiation. At log(λ) = 0.75 this community has split into two main smaller communities: the pink community contains constituent proteins of the RNA polymerase II mediator complex and the green community contains components of the closely related SAGA and TFIID complexes. At log(λ) = 1.6, this second community has split into the SAGA and TFIID complexes.

Multi-resolution community detection and characterisation is relevant both from the global viewpoint, where one can investigate the aggregate functional organisation of the proteome, and from the local perspective, where the community membership of particular proteins can be traced through changing resolution parameter. We thus now consider a protein-centred view of multi-resolution community detection. We consider, for an example protein, the properties of the communities to which it is assigned through changing resolution parameter, see Figure 5. The size of the communities, their mean similarity under the G and C measures, and the mean clustering coefficient are shown. The protein is a member of the ESCRT-I complex. (Figure S4 in Additional File 1 gives a further four examples.) Note the very robust properties of the communities in the A network over resolution parameter values of approximately 1 ≤ log(λ) ≤ 2.5, despite the tendency for them to be partitioned as λ increases. At these resolutions, the protein is in the same community as other members of the complex, as well as a few other very closely associated proteins. Beyond log(λ) = 2.5, the complex is broken up, as reflected in the drop in mean similarity values. The community present over 0.7 ≤ log(λ) ≤ 1.4 in the P network contains many proteins associated to the complex (in addition to the complex itself). Above the step observable at log(λ) = 1.4, only members of the complex are present. In Additional File 2, we give the names and brief functional descriptions of proteins that occur in some of the same communities for this example, and the four other examples given in Additional File 1. These five examples all show the following behaviour.

Use of topological properties to select functionally homogeneous communities

Almost all proteins are in functionally homogeneous communities at some value of the resolution parameter, and we therefore devise a method to swiftly identify these resolutions, especially if there is a dearth of functional information. We investigate whether any easily-calculated topological properties of the communities can act as indicators of functional homogeneity. Given a protein of interest we can then use such measures to quickly identify 'good' resolutions, without the need to assess functional homogeneity.

We find that the clustering coefficient is the most useful of the topological properties tested in the prediction of functional homogeneity for all three similarity measures and in both the A and P networks. The clustering coefficient of a network is a measure of the mean local clustering around nodes: A node has a high clustering coefficient, c, if its neighbours are also neighbours of each other [54, 55]. It is defined for each node as

(9)

where N_triangle is the number of triangles of which the node is a member, and N_triple is the number of connected triples of which the node is a member. (A connected triple is a single node with edges running to an unordered pair of other nodes.) Figure 6 shows for a) the A network and b) the P network the ROC curves for using the mean clustering coefficient of nodes in a community as a predictor of functional homogeneity for each of the three similarity measures in the A network. (See Methods for a description of the construction.)

There is some element of discretion for annotating A type interactions, i.e. deciding which pairs to list interactions between following experiments, with the principle competing models referred to as 'matrix' and 'spoke' [56]. This choice could cause artefactual topological features, so the extent to which we find particular topological features correlating with functional homogeneity could be sensitive to annotation choice. We are therefore encouraged that the same trends in predictive ability are evident in the P network, for which there is no such element of discretion.

As can be seen from Figure 5 and the figures in Additional File 1 Figure S4, clustering appears to be a good proxy for functional homogeneity when looking at individual proteins, and in the absence of much functional information could guide which resolution(s) should be targeted for investigation.