Abstract
Background
Formal classification of a large collection of protein structures aids the understanding of evolutionary relationships among them. Classifications involving manual steps, such as SCOP and CATH, face the challenge of increasing volume of available structures. Automatic methods such as FSSP or Dali Domain Dictionary, yield divergent classifications, for reasons not yet fully investigated. One possible reason is that the pairwise similarity scores used in automatic classification do not adequately reflect the judgments made in manual classification. Another possibility is the difference between manual and automatic classification procedures. We explore the degree to which these two factors might affect the final classification.
Results
We use DALI, SHEBA and VAST pairwise scores on the SCOP C class domains, to investigate a variety of hierarchical clustering procedures. The constructed dendrogram is cut in a variety of ways to produce a partition, which is compared to the SCOP fold classification.
Ward's method dendrograms led to partitions closest to the SCOP fold classification. Dendrogram or treecutting strategies fell into four categories according to the similarity of resulting partitions to the SCOP fold partition. Two strategies which optimize similarity to SCOP, gave an average of 72% true positives rate (TPR), at a 1% false positive rate. Cutting the largest size cluster at each step gave an average of 61% TPR which was one of the best strategies not making use of prior knowledge of SCOP. Cutting the longest branch at each step produced one of the worst strategies.
We also developed a method to detect irreducible differences between the best possible automatic partitions and SCOP, regardless of the cutting strategy. These differences are substantial. Visual examination of hardtoclassify proteins confirms our previous finding, that global structural similarity of domains is not the only criterion used in the SCOP classification.
Conclusion
Different clustering procedures give rise to different levels of agreement between automatic and manual protein classifications. None of the tested procedures completely eliminates the divergence between automatic and manual protein classifications. Achieving full agreement between these two approaches would apparently require additional information.
Background
This work investigates a large number of different methods for producing an automatic classification of structural domains of proteins based on allagainstall pairwise structural similarity scores. We produce candidate classifications and compare them to the human expertcurated classification SCOP [1]. Perhaps the earliest attempt at automatic classification of protein structures was by Holm and Sander [2], who produced the FSSP database, consisting primarily of a tree obtained by applying hierarchical clustering using the pairwise structural similarity scores among a sequence representative set [3]. Later, Dietmann et al [4] refined FSSP by introducing the notion of domain [5]. More recent work by Gewehr et al. [6] and others [710], focused on the problem of assigning new, unclassified domains into their correct preestablished fold classes. Other authors [1115] compared a set of pairwise similarity measures computed by various structure comparison methods to the SCOP (or other) preexisting classification, to investigate the causes of divergence between automatically determined pairwise structural similarity and expertcurated classification. In particular, the comparison between FSSP, SCOP and CATH, by Hadley and Jones [12], makes use of FSSP pairwise scores between protein chains which are present in all three databases, as FSSP does not use domain as the classification unit. Although these studies advanced the understanding of divergence between automatic similarity and expertcurated classification, improvement of automatic classification procedures was not their primary objective.
Our earlier work [15] pointed to structural variations within SCOP folds as a main cause of divergence between automatic and expertcurated classifications. Intrafold structural variation affects the measured structural similarity both within and between folds. The measured similarities are often not uniform among domains within a fold and the average similarity among domains within one fold can be different from that within another fold. We investigate different partitioning strategies in an attempt to accommodate such uneven similarity distributions. In particular, we investigate dendrogram cutting strategies as a potential means of isolating both tightly clustered, homogeneous folds and more heterogeneous ones, using a single procedure.
One source of potential confusion in the literature is the failure to explicitly consider the fundamental difference between mathematical properties of a similarity score matrix compared to a classification or partition. Mathematically, each can be represented as an M by M matrix where M is the number of protein domains. Unlike a pairwise similarity score matrix which can be any square matrix, a "partition matrix" is made of 0 or 1 elements indicating when two domains are in the same (1) or different (0) clusters of the partition. The partition matrix can always be transformed into a blockdiagonal form by sorting the rows and columns appropriately, reflecting a property related to the definition of a partition itself. A partition of a set of domains is comprised of nonoverlapping clusters, meaning that if A and B are in the same cluster, while B and C are also in the same cluster, then A and C are necessarily in that same cluster, a type of transitivity here called the "partition constraint". This property may be absent in a pairwise similarity matrix.
Building a classification by clustering based on pairwise similarity scores is essentially the same thing as transforming the pairwise similarity matrix into one which satisfies the partition constraint. Clearly, this process may force some domain pairs originally considered as dissimilar into the same cluster while other, similar pairs would be forced into different clusters. In many studies comparing the automatically determined pairwise similarities and a structural classification database [1115], the effect of the presence of the partition constraint was not explicitly considered. Some of the reported discrepancy between similarity scores and classifications may potentially result from the failure of the similarity scores to satisfy this constraint. In the present study, we eliminate this factor by first converting the pairwise similarity dataset into a partition (via clustering) and then comparing the two partitions.
The derivation of a partition from a set of pairwise similarities is not a trivial process, but involves several distinct steps of computations, requiring careful analysis. Here, particular attention will be given to hierarchical clustering methods and dendrogram cutting strategies. We obtain a partition from pairwise similarity scores by first transforming the similarity measure into a distance measure, applying various standard hierarchical clustering methods, obtaining a dendrogram, or binary cluster joining tree, and in the final step, applying a strategy which removes the root node and successively lower level nodes leaving behind a set of trees corresponding to the clusters in the partition. The pairwise similarity scores are provided by VAST [15,16], SHEBA [17] and DALI [18,19], three distinctive and efficient approaches for measuring structural similarity.
Results
Performance of hierarchical clustering as tree building methods
We generated pairwise similarity scores between all pairs of domains in the SCOP C class using three different structure alignment methods, VAST, SHEBA and DALI, and using two different similarity metrics for each method, altogether yielding 6 similarity matrices. We used four different hierarchical clustering methods (Single, Average, Complete linkage and Ward method [20]) to build dendrograms. The dendrograms were cut using seven SCOPindependent strategies and three SCOPdependent strategies. Finally, two different criteria for terminating the tree cutting process were investigated. We stopped the cutting either when a 1% false positive rate (FPR) was reached or when the number of clusters matched the number of folds in the SCOP C class, which is 94 for the current dataset. The true positive rate at 1% false positive rate (TPR_{01}) and the number of clusters in each partition, resulting from these various partitioning approaches, using the 1% FPR termination criterion, are reported in Tables 1 and 2.
Table 1. TPR values at 1% FPR (TPR_{01}) for Hierarchical Cluster Methods and Tree Cutting Strategies.
Table 2. Number of Clusters for Hierarchical Cluster Methods and Tree Cutting Strategies at 1%FPR.
Trees generated by Ward's method result in better performance. The TPR_{01 }values obtained for Ward's clustering method are, with one exception, always higher than for the Average, Complete or Single linkage method. The simplest, level cut strategy with Ward's method achieves an average TPR_{01 }of 50% across all six structure comparison scores. Complete or Average linkage clustering constitute the next best alternatives to Ward's method, with average TPR_{01 }across all similarity scores and structure comparison methods, of 32% and 29% respectively. Single linkage clustering shows uniformly the most divergence from SCOP. Not only does Ward's method achieve the highest average TPR_{01 }value, but the values vary less across different similarity scores and structure comparison methods, than do those for Complete, Average and Single Linkage, suggesting that Ward's is more satisfactory for this application.
With other tree cutting strategies the trend among various clustering methods is the same. For the largest size cut strategy, Ward's method gives an average TPR_{01 }value of 61%, at least 10% better than for the three other clustering methods. Again, Complete and Average Linkage give similar average TPR_{01 }values, while Single Linkage is very low in comparison. In view of the clearly superior results with Ward's method, we investigate tree cutting strategies using that method alone.
Performance of tree cutting strategies
In the following, the performance of treecutting strategies is analyzed. The 10 tested treecutting strategies fall into four groups according to their TPR_{01 }values (Figure 1 and Table 1), and are represented by the longest branch cut (1), level cut (2), largest size cut (3) and mutual information cut (4) strategies. Figure 1 illustrates the comparative performance among these representative tree cutting strategies, for the VAST number of matched residues similarity score.
Figure 1. Receiver Operating Characteristic (ROC) curves for various tree cutting strategies, using the Ward's method clustering based on the VAST Nres similarity score. The true positive rate (TPR, Y axis) is plotted against the false positive rate (FPR, X axis). Black curves correspond to largest size cut (solid line), longest branch cut (dotted line) and level cut (dashdot line) strategies, which use only of tree topology parameters. The red dotted curve corresponds to the MI cut strategy, which uses SCOP classification information directly, so is viewed only as an approximate upper bound to the SCOPindependent strategies. The Yaxis intercept for this curve is 0.45 and close to 0.0 for all SCOPindependent strategies. Only the portion of the curve for FPR values below 0.45 is shown. Between this FPR value and 0.63 relative positions of the curves do not vary. Above 0.63 FPR value all curve superpose to the MI cut curve.
Among the four SCOPindependent tree cutting strategies shown, "largest size cut" achieves the best performance with its ROC curve dominating all other SCOPindependent ROC curves. Moreover, "largest size cut" ROC curve is the closest to that for Mutual Information cut, a strategy which attempts to maximize the agreement of the partition with SCOP at each step. Similar results (data not shown) were seen for SHEBA and DALI similarity scores.
The longest branch cut group (Table 1, Group 1) contains the strategies of lowest performance, e.g. tree skewness strategy, as indicated by their respective mean TPR_{01 }values. Longest branch cut strategy results in one of the lowest agreements with SCOP, in particular for FPR values below 30% (Figure 1). Above this FPR level, its performance becomes more acceptable relatively to the upper bound given by the mutual information cut. This strategy seems to make more sense when intercluster distance is high, i.e. for highest nodes of the tree. But it behaves much worse than others when the intercluster distance is reduced, as its ROC curve is closer to the main diagonal than the ROC for other strategies. In the lowest part of the tree, intercluster distances vary only slightly so that the longest branch criterion for choosing among the many possibilities, the next subtree to be cut might not be discriminative enough.
Level cut and maximum entropy cut strategies form another group (Table 1, Group 2) characterized by middling performance. Their grouping is explained by the fact that these strategies perform comparably at low FPR values and their average TPR_{01 }values are comparable. Higher variability among the TPR values for maximum entropy cut strategy is noted, however. Indeed, for VAST and SHEBA, level cut strategy is better than maximum entropy cut, as a rule, with DALI number of matched residues (Nres) score the only exception. The highest TPR_{01 }for the maximum entropy cut strategy (58%) is obtained with the Nres metric, while the lowest value (40%) is obtained by VAST Pcli.
The third, largest size cut strategy group (Table 1, Group 3) also includes the tree completeness cut, and highest tree cut strategies and corresponds to SCOPindependent strategies of highest performance, with average TPR_{01 }values ranging from 58 to 61%. The tree height and tree completeness are closely related to the number of leaves, i.e. the size of the tree. These results indicate that of these three topological properties, the size is a better measure, for the purpose of creating partitions that agree with SCOP. This is true for partitions of small size (the large FPR values), but also when the number of clusters is large (the small FPR values).
The fourth group (Table 1, Group 4) is made up of mutual information cut and Best TPR/FPR Ratio cut strategies, all of which make direct use of the SCOP fold partition and thus are not useful for an independent classification. The poor performance of the best TPR cut strategy (Table 1, Group 1) indicates that strategies which ignore the false positives will not generate SCOPlike partitions.
The size of each partition is reported in Table 2, and shows a large variation in the number of clusters depending on the hierarchical clustering method. Single linkage always leads to partitions with a high number of clusters, in fact close to the total number of domains M, when FPR is kept low at 1%.
Comparison of partitions
The relative organization of automatically generated partitions can be understood by first finding a distance measure between two partitions and then displaying the partitions in a distance preserving graph. The distance between partitions (Eq. 8), computed as the total number of disagreements about whether or not a pair of domains is in the same cluster, is reported in the Table 3. In this exercise, partitions with exactly 94 clusters were compared, including the SCOP fold partition.
Table 3. Distance between partitions*, (Δdistance in 1000s)
Average distances (Table 4) among partitions within a given comparison method (either VAST, SHEBA or DALI), are uniformly lower than distances between these latter and SCOP, by almost a factor two. Partitions from VAST tend to be slightly more heterogeneous than those from DALI or SHEBA. On average, automatic partitions from the three comparison methods are similarly distant from SCOP. Similarly, for a given cut strategy, the average distance among its partitions is consistently smaller than the distances from those automated partitions to SCOP (Table 4). Further, largest size cut partitions are half as distant among themselves as are those of level cut strategy. The pattern of distances confirm that the largest size cut partitions are much closer to the SCOP fold partition than level cut partitions, as also seen using the ROC curves.
Table 4. Δdistance*, in 1000s, between partitions within each method and minimum Δdistance from each method to SCOP
Table 5 shows that every automatic partition is much closer to other automatic partitions, than to the SCOP fold partition. For example, the DALI Zscore levelcut partition (DZL) is at most 56,000 units from all other automated partitions, but 91,000 units from SCOP. The closest automated partition to SCOP is DALI Zscore, largest size cut (DZS), with a distance of 63,000 units, which in turn is no farther than 49,000 units from all other automated partitions. DZL is farthest from both SCOP and from at least one other automatic partition. Conversely, DALI Nres with Largest sizecut (DNS) is perhaps most representative of automated methods as it minimizes the maximum distance to other such methods.
Table 5. Distance between Automatic Partitions and SCOP fold Partition
A convenient, intuitive representation of the organization of these partitions is obtained using multidimensional scaling, a technique which embeds the 13 partitions into a lowdimension Euclidean space so that the pairwise distances are approximately preserved. Figure 2 represents the automatic and SCOP fold partitions in a 2dimensional space. Automatic partitions are wellseparated from the SCOP fold partition, which appears isolated. Largest sizecut partitions are generally closer to SCOP than are levelcut partitions, and are also less spread.
Figure 2. Classical multidimensional scaling (MDS) plot of various automatically generated partitions and the SCOP fold partition. The MDS plot approximately preserves distances, Δ, between partitions (see Methods). The X axis represents the projection of the location of each partition onto the first eigenvector (first principal component), while the Y axis is a projection onto the second eigenvector (second principal component). The X and Y axes of the plot are scaled to represent the distance, Δ, divided by 1000. Automatic partitions are obtained with Ward's method clustering, based on two different similarity score for each of the three methods VAST, SHEBA and DALI. Each partition is designated by an uppercase letter: A, VAST Pcli Level cut; B, VAST Nres Level cut; C, DALI Zscore Level cut; D, DALI Nres Level cut; E, SHEBA Zscore Level cut; F, SHEBA Nres Level cut; G, VAST Pcli Largest Size cut; H, VAST Nres Largest Size cut; I, DALI Zscore Largest Size cut; J, DALI Nres Largest Size cut; K, SHEBA Zscore Largest Size cut; L, SHEBA Nres Largest Size cut, and S, the expertcurated partition SCOP. Automatic partitions resulting from the same treecutting strategy, are grouped together within the same dotted area.
Comparison of cluster size distributions
In addition to the number of clusters in a partition, the distribution of cluster sizes may be of interest in selecting an appropriate classification. Figure 3 shows the cluster size distribution for six automatically generated partitions and for SCOP. For comparison, a partition of the same number of domains randomly assigned to 94 clusters with equal probability is shown. This distribution was approximated by a Poisson distribution with a mean value of 1330/94 ≅ 14.15.
Figure 3. Box plots of the cluster size distributions for six automatically generated partitions with 94 clusters and the SCOP partition of the C class. Each partition is associated with a box plot. A box plot summarizes the following statistics: median (line within the box), upper and lower quartiles (the upper and lower hinge of the box respectively), minimum and maximum data values (the ends of the vertical dash lines), and outliers (circles). For comparison, a "random" partition would have cluster size following approximately a Poisson distribution with intensity parameter equal to M/94 ~ = 14 domains per cluster. The partitions are labelled as explained in Figure 2. The SCOP partition and Random partitions are identified. All distributions necessarily have the same mean value, M/94, but show differing medians, interquartile ranges, tail lengths and maximal values. Both the SCOP and Random distributions show minimum values substantially above the automated partitions. The horizontal red dash line indicates the median value of SCOP fold partition.
We observe first that all partitions, including SCOP, have lower median cluster size and greater spread of size than for random. There is evidently sufficient signal strength within the similarity score matrix to influence the size distribution. Second, the 75th percentiles for automated methods tend to be larger than for SCOP, while the SCOP distribution shows larger positive skewness, with a greater number of unusually large clusters. Third, there is some uniformity in the size distribution within tree cutting strategies, with largest size cut showing somewhat higher 75th percentiles and less skewness than do the level cut distributions. Level cut distributions are closer the SCOP distribution in terms of median, 75th percentile size and the larger number of outliers, than the largest size cut distributions.
The largest SCOP fold (c.1), is in fact, split by all strategies and methods as none include a cluster with 182 domains. Largest size cut strategy intentionally eliminates outliers of large size, thereby creating more clusters of intermediate sizes, with greater spread (interquartile range) than for level cut. The level cut generates a few large size outliers but the clusters are smaller, typically. This is consistent with the observation made earlier that the level cut behaves like the maximum entropy cut at small FPR ranges. Using partitions of 94 clusters, we find FPR values higher than 1% for level cut strategy but lower than 1% for largest size cut strategy (Table 5), indicating that a tradeoff must be made in matching the size distribution as well as maximizing the TPR in relationship to the SCOP partition.
Identification of dispersed folds
The spanning cluster of a SCOP fold is the smallest cluster in the dendrogram which spans or includes all domains in that fold. The excess span of that fold (see Methods) are the domains from other folds that are included in its spanning cluster. A homogeneous cluster is a cluster which includes only domains from a single SCOP fold. The size of the excess span and the size of the largest homogeneous cluster are given in Table 6 for each fold in SCOP C Class for three different dendrograms. These two measures allow comparison of each dendrogram to SCOP on a fold by fold basis, and can highlight regions of agreement or disagreement between the two systems. When the excess span is zero and the largest homogenous cluster is 100% of the fold size, the dendrogram and SCOP are in perfect agreement for that fold. A hypothetical treecutting strategy could potentially isolate this particular subtree to form a cluster exactly matching that fold.
Table 6. Excess Span and Largest Homogeneous Cluster Size for SCOP C Class folds.
The first 38 folds reported in Table 6 are in perfect agreement with dendrograms of all three structure comparison methods VAST, SHEBA and DALI. Together they comprise 187 of 1330 domains in the Cclass. The next 23 folds (Table 6, rows 39–61, comprising 227 domains) agree perfectly with the dendrogram of at least one structure comparison method. Thus 61 out of 94 SCOP folds in the C class are consistent with the dendrogram built from at least one pairwise structural similarity measure.
Structure comparison methods differ in the consistency of their associated dendrogram with the SCOP folds. Dendrogram derived from each similarity score method disagree with SCOP for various folds. Fortyeight folds (comprising 1080 domains) disagree for the VAST tree, 41(1026 domains) for the SHEBA tree, and 43 (991 domains) for the DALI tree, according to this criterion.
We consider a fold to be highly dispersed if it disagreed with trees of all three structure comparison methods. There are 33 such folds (comprising 916 domains) and they are reported in Table 6, row 62–94. None of these 33 folds could be obtained as a homogeneous cluster, thus each contributes to the loss of agreement between any automatic partition and the SCOP fold partition. For these 33 folds, the intersection of their excess span was computed, and reported in Table 4, column labelled INT. INT is a measure of the disagreement with SCOP that remains even if all three dendrograms were combined. Thirtyone out of 33 folds give rise to a positive INT, meaning that same domains contributed to their dispersion within trees. These 31 folds therefore contribute to the remaining distance between the automatic and expertcurated partitions, regardless of tree cutting strategy or structure comparison method used.
Dispersion caused by low structural similarity within folds
To examine such dispersed folds in detail, we select two examples. Figure 4 schematically represents the situation of fold c.58 within each dendrogram. Instead of forming one homogeneous cluster including all domains of the fold, c.58 consists of mainly two homogeneous clusters set far apart in the tree. Between these two homogeneous clusters, domains from other C class folds intervene, in particular from c.78. Thus, one homogeneous cluster of fold c.58 is considered more similar to domains from fold c.78 than to other domains from fold c.58, and this situation pertains to all three structure comparison methods.
Figure 4. Schematic of dispersion of fold c.58. Nodes labelled c.58 represent two of the homogeneous clusters of fold c.58. Node labelled c.78 is one homogeneous cluster of fold c.78. Relative join positions reflected from the complete Ward's method dendrograms, based on clustering the Nres similarity scores.
Figure 5a, b and 5c presents the matrices of pairwise distances between all domains within folds c.58, for VAST, SHEBA and DALI, respectively. Coding the distances by color makes obvious why fold c.58 is split into mainly two homogeneous clusters by all three dendrograms. Figure 5a, b and 5c also include domain d1b74a1 from fold c.78, which is highly similar to domains of fold c.58. Indeed, many of the pairwise distances involving d1b74a1 (coded yellow), are small enough that it would fit into any cluster of fold c.58. This pattern of high distance between homogeneous clusters of c.58 and low distance of a fold c.78 domain to members of c.58 explains why no clustering method and no treecutting strategy is likely to perfectly identify fold c.58.
Figure 5. Matrices of pairwise distances among domains of folds c.58 and c.78. Pairwise distances between domains were obtained from the pairwise similarity scores number of matched residues, according to Eq(1) and symmetrized as specified in Methods. The yellow color corresponds to pairwise distances smaller than141, 119 and 88 for VAST matrix (a), SHEBA matrix (b) and DALI matrix(c) respectively. The white color corresponds to pairwise distances higher than 141, 129 and 88 for (a), (b) and (c) respectively.
This pattern of structural similarity measures is further analyzed based on structural superposition of domains. Figure 6 (a) and 6 (b) superpose domain pairs from the same homogeneous cluster, and show the degree of structural similarity typical within cluster. The structural superposition of the domains from different homogeneous clusters of c.58 is shown in Figure 6 (c). The low similarity observed here is mainly due to the difference of the N terminal features, with the 3 layer a/b/a feature typical of this fold present in both domains. Figure 6 (d) is the structural superposition of a domain from c.58 with one from c.78. Again, the strands and helices making the 3 layer a/b/a feature are well aligned.
Figure 6. Superposition of c.58 and c.78 domain pairs, by DALI alignment. Panel (a): one intracluster, cluster 1, pair of c.58, d1edza2 (146 residues, red) and d1b0aa2 (121 residues, blue), 118 residues aligned, RMSD = 2.7Å, and Zscore = 14.4. Panel (b): another intracluster, cluster 2, pair of c.58, d1c1da2 (148 residues, yellow) and d1leha2 (134 residues, magenta), 130 residues aligned, RMSD = 1.4Å and Zscore = 21.3. Panel (c): intercluster pair of c.58, between cluster 1 and cluster 2, d1b0aa2 (blue) and d1c1da2 (yellow), 81 residues aligned, RMSD = 4.7Å and Zscore = 4.5. Panel (d): interfold domain pair between d1b74a1 (c.78, 105 residues, green) and d1b0aa2 (c.58, 121 residue, blue), 75 residues aligned, RMSD = 3.3Å, Zscore = 6.2. Both c.58 and c.78 folds are characterized by having a core of 3 layers a/b/a with a parallel betasheet of 4 strands, ordered 2134. Under hierarchical clustering, fold c.58 is comprised of two homogeneous clusters (Figure 4).
In Figure 6, the structural superposition (d) of domains from different folds produces even a better alignment than the superposition (c) of domains from two different clusters within the same c.58 fold. The disagreement between automatic and expertcurated classification, arises directly from the low similarity within fold c.58, and the substantial similarity to domains of another fold. It becomes impossible to merge the two distinct homogeneous clusters of fold c.58 without including domains from c.78 as well.
As another example, we select fold c.1, one of the largest in the dataset, and a highly dispersed fold. The intersection of spanning clusters for this fold includes three domains from the single fold c.6, Cellulases, which are partial barrels. Fold c.6 is found to be easily identifiable by all three similarity methods (Table 6, row 22). Figure 7 shows two TIM barrel fold structures and one Cellulases structure. TIM barrel structures are easily recognizable and are not likely to be confounded with any other folds by the human expert. The typical TIM barrel structure d1clxa_ (Figure 7 a) appears to be more similar to the Cellulases structure d1dysa_ (Figure 7 b), than to another TIM barrel, d1a4ma_ (Figure 7 c), by all three structure comparison methods. TIM barrel structure (c) has 8 strands although two of them are much longer and two are much shorter than the rest making the barrel somewhat distorted, compared to the typical TIM barrel structure (a). The lower similarity between a typical TIM barrel domain and another member of that fold, compared to the similarity of a c.6 domain to the typical TIM barrel, means that automated clustering methods cannot separate these two folds perfectly. The structural distinctions between the two folds are evidently too subtle to be detected by these structure comparison methods.
Figure 7. Structures from TIM barrel fold c.1 and its variant fold c.6. Panel (a): d1clxa_, a typical TIM Barrel domain belonging to family c.1.8.3 with 345 residues. Panel (b): d1dysa_ c.6.1.1 with 345 residues. Panel (c): d1a4ma_, also a TIM Barrel, belonging to family c.1.9.1 with 349 residues. These domains belong to three different homogeneous clusters, by all three programs. The pairwise distances by VAST, SHEBA and DALI, between structures (a) and (c) are (334*, 440, 318), between (b) and (a) are (236, 350, 256), and between (b) and (c) are (420, 500, 344). The pairwise number of matched residues (Nres) by VAST, SHEBA and DALI between structures (a) and (c) are (180*, 127, 188), between (a) and (b) are (227, 170, 217) and between (b) and (c) are (137, 97, 175). *This corresponds to the best alignment of the two non symmetric pairwise alignments produced by VAST.
Discussion
There are two approaches for comparing results from automatic structure comparison methods to an expertcurated reference classification such as SCOP. One can either directly compare the pairwise structural similarity measures to a similar measure derived from the reference partition, or produce a partition from the pairwise structural similarity score and compare it to the reference partition. We adopted the latter approach in this study, as we expected that enforcing partitioning constraints would introduce a new element that is not present in the pairwise similarity measures alone.
We explored a variety of methods for obtaining automatic partitions from pairwise structural similarity measures computed by VAST, SHEBA and DALI. Specifically, four different hierarchical clustering methods were used to construct dendrograms or binary trees, which were then cut into subtrees by ten different treecutting strategies, to produce partitions. The results show that the combination of Ward's method with the largest size cut strategy has best agreement with SCOP among all combinations of clustering methods and treecutting strategies explored so far.
Trees generated by Ward's method result in partitions agreeing better with SCOP folds, than those generated by Single, Complete or Average linkage, regardless of the tree cutting strategy applied. Clusters formed by Ward's method tend to be highly concentrated around a mean, as they are formed so that the variance within the cluster is minimized. This suggests that SCOP folds are constructed based on the cohesiveness of the group as a whole rather than on similarity of individual pairs. This may be a more appropriate view of folds than as complete subgraphs where every structure is related to every other structure as suggested by Complete linkage clustering. Although Average linkage produces clusters which are organized around a constructed mean, it does not minimize the variance around this mean, in contrast to Ward's method. Its lower performance relative to Ward's method strengthens the view of folds as structurally cohesive groups of domains. The mediocre results obtained by Single linkage emphasize even more the suitability of the notion of structural cohesiveness in modelling folds. Indeed Single linkage folds can be seen as minimum spanning trees, which do not necessarily require structural cohesiveness.
Largest size tree cutting strategy and those associated with it, such as highest tree cut and tree completeness cut, achieve much better agreement with SCOP, than other SCOP independent strategies. In particular they outperform the "standard" level cut strategy. This is an unexpected result and suggests that the fold size was perhaps an implicit criterion, among others, in the formation of SCOP folds and that large groups of protein structures were split even when their internal similarity was higher than the internal similarity of smaller sized group of proteins. The small spread of the sizes of SCOP folds shown in Figure 3 supports this notion. On the other hand, a highly structurally distinctive fold such as TIM barrel (c.1) which is highly populated relative to the majority of folds and easily, visually recognizable, form an exceptionally large cluster.
The best combination of hierarchical clustering and SCOP independent tree cutting strategy e.g. Ward's method with either largest size, tree completeness or highest tree cut, resulted in an average TPR_{01 }values (61%, 58% and 59% respectively) that are only slightly above the average agreement of 57% achieved by comparing pairwise similarities directly, without clustering, to SCOP folds partition, using the same set of SCOP C class domains. This suggests the idea that the persistent discrepancy between automatically determined similarity and SCOP is most likely not due to the partitioning constraint, even though all possible partitioning strategies, such as most recently developed clustering techniques in [21,22] for example, have not been examined. Figure 2 provides additional evidence supporting this view. It shows that automatic partitions produced by introducing partitioning constraints into pairwise structural similarity measures, cluster together in the space of partitions, far away from the expertcurated classification SCOP implying that automatic structure comparison methods agree well with each other, and that there might be irreducible differences between them and SCOP.
Automatic methods VAST, SHEBA and DALI obtained 71%–76% TPR_{01 }using the combination of Ward's method and the MI cut strategy which maximizes the agreement with SCOP. This is higher than the maximum 65% and 67% obtained by the largest size cut and the direct pairwise comparison scheme, respectively, on the same set of the SCOP C class protein domains. As the MI cut strategy proceeds with prior knowledge of SCOP folds information, this range of agreement is close to the maximum obtainable between automatic classification and manually curated SCOP dataset, regardless of the partitioning procedure. The inherent level of disagreement can be seen from the number of SCOP folds for which the spanning cluster exceeds the size of the fold, in Table 6.
Examination of two such folds corroborates the findings of our previous paper [15] that some SCOP folds includes structural variation causing measured similarity between members of the same fold to fall below that between members of related, but different folds. Indeed, the distance matrix in Figure 5 shows that fold c.58 displays this pattern, due to low similarity among domains from its distinct homogeneous clusters, and a high similarity with domains from a different fold c.78. The superpositions in Figure 6 (a), (b) and 6 (c), indicate that in all cases, be it intra or intercluster, all three structure comparison methods aligned the structural feature defining the fold, and yet intracluster structural similarity is higher than intercluster structural similarity. Regarding the TIM barrel fold, a detailed study by Nagano et al. [23] also characterized the fold as highly diverse, although using the CATH [24] database. Thus, structural variation within folds, shown in [15] as the main cause of divergence between automatic and expertcurated classifications, does not disappear even after satisfying the partitioning constraint. This strengthens our previous findings and emphasizes the importance of properly handling the structural variation in order to reduce the gap between automatic and expertcurated partitions. We also have underlined some relationships between distinct folds. The domains which correspond to the excess span of a dispersed fold could be seen as evidence suggestive of a evolutionary relationship among folds, as discussed by Lupas et al [25]. Domains from fold c.6, which are excess span of TIM barrel fold, are actually variants of TIM barrels.
Finally, in the light of this analysis, we think that future improvements to automatic protein structures classification would likely come by explicitly identifying common structural cores. Pairwise similarity scores alone appear to be limited in that regard, so techniques involving multiple structural alignment will likely be needed.
Conclusion
The level of agreement between manual and automatic classifications varies with clustering methods and tree partitioning strategies. However, the best agreement reaches similar upper bound than when structural pairwise similarity is compared directly to the manual classification. Therefore divergence between automatic and manual classifications is not eliminated by the introduction of partitioning constraints.
Our observations are based on SCOP C class, but are likely valid for other classes as well, as C class domains contain both types of secondary structure elements, alpha helix and beta strand, and are the most difficult to classify due to higher confusion among them compared to the all alpha or all beta classes [15]. Our exploration of potential classification procedures of proteins based on structural similarity is complementary to the analysis done in our previous paper [15]. The modular structure of the proteins has been accounted for by using a database of structural domains defined in SCOP although this domain parsing may be at variance with other domain parsing such as that in CATH [24], or those based on amino acid sequences in CDD[26], ProDOM [27] or Pfam [28]. We did not address this issue here. Notwithstanding efforts made in this present work, there is still a discrepancy between the results of the automatic structure comparison methods and the SCOP classification. Based on the maximum attainable TPR of about 76%, roughly one quarter of the C class domains are not classified by SCOP according to measured global structural similarity. As we previously suggested, the main reason is that global structure similarity cannot entirely account for the characteristic local structural features on which the SCOP classification is based. Future research should be aimed at finding algorithms able to automatically extract such evolutionarily conserved, common local structural features of domains. As observed by Chothia and Lesk [29], in proteins having low sequence identity and the same biochemical function, only about half of secondary structures are conserved. The remaining challenge is to identify the conserved half!
Methods
We consider a set E of M = 1330 domains in the C class selected from the set of SCOP domains with less than 40% pairwise sequence identity in ASTRAL [30] version 1.63. Their associated N = 94 folds constitute a partition, i.e. a collection of disjoint folds or subsets of E whose union is E.
Our objective is to compare this fold partition, an expertcurated classification, with a partition of the same set E of domains, obtained automatically as follows.
First, all M^{2 }= 1,768,900 pairs of domains drawn from this set E were compared by VAST, SHEBA and DALI. The calculations were made using the highperformance computational capabilities of the Biowulf PC/Linux cluster at the National Institutes of Health, Bethesda, MD [31].
Data resulting from the computation were used to produce matrices of similarity scores [S_{qt}], where S_{qt }is either the Nres (Number of Matched Residues) between domains q and t, or the Zscore for SHEBA and DALI and the Pcli score for VAST.
Second, these similarity matrices were transformed into distance matrices and symmetrized by replacing the upper diagonal values with the lower diagonal values. There exists many ways of obtaining a pseudodistance measure from a similarity measure [32], even though the triangular inequality is not always guaranteed. We choose a simple approach which guarantees that selfdistance is zero and all distances are nonnegative. With similarity defined as the number of matched residues, the distance satisfies the triangular inequality. For domains q and t within E, and similarity measure S_{qt}, the distance D_{qt }is defined as follow:
These matrices of pairwise distances were then entered into hierarchical clustering algorithms to obtain dendrograms or binary trees. We considered Matlab implementation [33] of four representative hierarchical clustering methods [20]: Single, Average, Complete linkage and Ward's method.
Tree cutting algorithm
A dendrogram resulting from a hierarchical clustering is a binary tree, where each node is associated with an intercluster distance defined by the joining distance between its left and right children. The intercluster distance for the leaves of the tree is defined to be zero. To obtain a partition from a dendrogram, one ordinarily chooses a level, i.e. an intercluster distance value, falling within the dendrogram, then removes or "cuts" all the joinnodes above the chosen level, leaving behind a set of trees whose roots or joinnodes fall below that specified level. Here, we generalize this approach and find partitions made up of a set of subtrees which join with root nodes at various levels, by introducing a panel of tree cutting criteria making use of a variety of tree characteristics, in addition to the intercluster distance. Our recursive algorithm starts with the original tree. When its root is cut or deleted, two trees remain, representing the two clusters in the growing current partition. Depending on the number of clusters desired and various topological features of each of the remaining trees, one tree is selected and its root is cut, leaving behind two smaller trees, each again representing a cluster in the refined partition. The procedure is terminated when the desired stopping criterion is reached.
More formally the tree cutting algorithm starts from the initial partition consisting of one cluster containing all elements of E. Given a partition of size K > 1, a partition of size K + 1 is obtained from K by splitting one of the clusters K_{i }into 2 distinct clusters. The splitting of a cluster is represented by removing the root node of its associated tree.
At the initial step of this process, there is no choice, but cutting the root node of the tree. Further partitioning involves a choice of which cluster to split. A tree cutting strategy determines this choice. We present below several strategies which have been developed here and used in this analysis. We consider three different ways of stopping the tree partitioning process. First, the partitioning process stops when no further partition can be obtained, corresponding to the situation where all clusters of the partition are singletons. Second, the partitioning process stops when a partition with a given number of clusters is obtained. We will be mainly interested by partitions of size 94 clusters corresponding to the number of SCOP C class folds in the dataset. Third the partitioning process stops based on reaching < = 1% FPR (see below) when comparing the automatic partition with SCOP.
Tree cutting strategies independent of SCOP
We define seven SCOP independent tree cutting strategies. Such a strategy determines the next subtree to cut, without using a prior knowledge of SCOP fold partition. They all have been implemented for this analysis.
Level cut
The strategy proceeds by descending intercluster distance, starting at the level of the root node. At a given step of the tree cutting strategy, resulting trees are available candidates for further cutting. The level of the next cut is determined by the tree whose node has the highest intercluster distance among all candidates. An illustration of the level cut is given in Figure 8.
Figure 8. Schematic of cutting strategies and spanning clusters. Leaves of the binary tree, e.g. domains, are numbered 1 to 12. In the figure, it is assumed that four clusters (1,2,3), (7, 8, 9) and (12) are from one particular fold and the three clusters (4, 5, 6), (10) and (11) contain all 5 domains of another fold. Horizontal line A – A represents the level cut, which produces a partition of 7 clusters but which splits the cluster (1,2,3) into two clusters, (1,2) and (3). Oblique line B – B is an allowable cut which produces a partition of 6 clusters and does not split (1,2) from (3). The node n1 represents the spanning cluster of the fold having domains 4, 5, 6, 10 and 11. The span, or size, of this spanning cluster is nine, and its excess span, or number of included domains not in this fold, is four. The node labelled n2, spanning three domains, namely 4, 5, and 6, represents the largest homogeneous cluster of this fold. For this fold, the relative size of the largest homogeneous cluster is 3/5, e.g. the size of the largest homogeneous cluster divided by the size of the associated fold.
Largest size cut
At a given step of the tree cutting strategy, resulting trees are available candidates for cutting. The candidate tree which contains the largest number of domains or leaves, is cut. This strategy counts the exact number of domains in the subtree evaluated for cutting.
Highest tree cut
The topological height h(K_{i}) of the node K_{i }is the number of intermediate nodes along the longest path from K_{i }to the leaves. The leaves are of topological height zero. For each cluster K_{i}, the topological height of its associated subtree is computed. The tree with the topologically highest root node is cut.
Tree completeness cut
Given the height h(K_{i}) of the subtree with root node K_{i}, the number of nodes in this subtree if it were complete would be 2^{h(Ki)+1}1. If the subtree is not complete, its actual number of nodes, including the root and the leaves, is less than 2^{h(Ki)+1}1. Tree completeness is defined as the ratio of the actual number of nodes of the subtree, and the number of nodes if it were complete. The smaller the ratio, the less complete is the tree. The tree with the smallest ratio (least complete) is cut.
Tree skewness cut
Given a tree, its skewness is defined as the ratio of the number of nodes in its left and right children. The ratio is defined such that the number of nodes of the smallest child, left or right, is divided by the number of nodes of the largest child. The tree of greatest skewness is chosen for cutting.
Longest branch cut
The branch length of a tree is the difference between the intercluster distance of the root node and the smaller of the intercluster distances of its two children. The tree with the longest branch value is cut.
Maximum entropy cut
The entropy H(K) of an entire partition K is defined by[34]:
where , the probability of the cluster K_{i}, is ratio of the number of domains within K_{i }to the total number M of domains being partitioned. The tree which would result in the greatest increase in entropy value for the current partition is chosen for cutting.
Tree cutting strategies which refer to SCOP
We define here three tree cutting strategies which use the prior knowledge of SCOP fold partition to determine the next subtree to cut:
Mutual Information (MI) cut
The mutual information I(K:F) between a partition K resulting from the cutting of a binary tree and the SCOP fold partition F, is defined as follow:
with P_{i }as defined for entropy, and is ratio of the number of domains in common to cluster K_{i }and fold F_{j }to M. Cutting a particular tree results in a new partition, as well as a new value of I(K:F). The tree which would result in the greatest increase of the mutual information is cut.
Best TPR/FPR Ratio cut
TPR and FPR values (defined below) refer to the SCOP fold partition. The tree whose cutting would results in the greatest TPR/FPR ratio, is cut.
Best TPR cut
The tree which would result in the greatest TPR value is cut.
ROC curves for comparing partitions
The ROC curves are constructed using the definition of true and false positive rates introduced in [15], but with the notion of membership in a cluster replacing the notion of similarity cutoff value. The following definitions reflect this variation.
The true positive rate between two partitions is defined as
where TPR_{j }is defined by:
where n_{jk }is the number of domains common to cluster k and SCOP fold j, and n_{j }is the number of domains in fold j.
The false positive rate between partitions is defined as follow:
where
when i and j referring to SCOP folds i and j.
Distance between partitions
Given two partitions K and P, we define the distance between these two partitions as follow:
where K and P are binary incidence matrices of the two partitions whose elements take on the value 1 when the domains i and j are in the same cluster, and 0 otherwise.
Multidimensional scaling (MDS) plot
We use classical multidimensional scaling to plot a set of points in Euclidean space, such that each point represents a partition and the Euclidean distance between two points in the plot approximates the value of the Δdistance between partitions they represent. Thirteen partitions are represented in the plot, including the expertcurated fold partition SCOP. Automatic partitions were chosen to have as many clusters as folds in the C class, by requiring the partitioning algorithm to stop when a partition of size 94 clusters is reached. The computation of pairwise Δdistance between the 13 partitions resulted in a 13 by 13 matrix (see Table 3) which is then used as input to the classical MDS procedure in Matlab [35]. A two dimensional plot of the 13 partitions is obtained by projecting the points onto the two X and Y axes where the X axis is the eigenvector corresponding to the largest positive eigenvalue, and the Y axis the eigenvector corresponding to the second largest eigenvalue.
Measure of fold dispersion
We define a measure of the dispersion of a SCOP fold within a binary tree, based on two metrics: the size of the excess span of the fold, and the size of the largest homogeneous cluster of the fold. The spanning cluster of a SCOP fold is the smallest cluster in the tree including all its domains. Proceeding upward from the leaves (domains) of that fold and stopping at the first (lowest) node common to all its leaves, we find the spanning cluster of the fold. The size of this spanning cluster, i.e. the number of leaves or domains it includes, is always greater or equal to the size of the associated fold. The set theoretic difference between the spanning cluster and the associated fold, is called excess span. The intersection of spanning clusters from distinct trees formed using distinct methods, identifies domains which might reasonably be considered as similar to domains in that fold. The size of the intersection less the members of that fold is a measure of how much larger the fold should be if it were to include all reasonably similar domains.
A homogeneous cluster is one containing domains of only one fold. Folds may be comprised of more than one homogeneous cluster. The size of the largest homogeneous cluster as a fraction of the size of a fold measures the fraction of the fold which can be easily recognized by the similarity measure and clustering algorithm. Figure 8 illustrates the notion of homogeneous and spanning clusters for a schematic binary tree.
Authors' contributions
VS, PM – execution of pairwise comparisons using VAST on Biowulf computer, development of partition algorithms, statistical analysis
JFG, JG – development of VAST program, interpretation of results
CHT, BL – development of SHEBA program, and pairwise comparisons, 3D structure visualization, interpretation of results
All authors participated in the formulation of the problem and in devising specific solutions to the problem. All authors have read and approved the final version of the manuscript.
Acknowledgements
Molecular graphics images were produced using the UCSF Chimera package [36] from the Resource for Biocomputing, Visualization, and Informatics at the University of California, San Francisco (supported by NIH P41 RR01081).
This research was supported in part by the Intramural Research Program of the NIH, Center for Information Technology and the National Cancer Institute.
References

Murzin AG, Brenner SE, Hubbard T, Chothia C: SCOP: a structural classification of proteins database for investigation of sequences and structures.
Journal of Molecular Biology 1995, 247:536540. PubMed Abstract  Publisher Full Text

Holm L, Sander C: The FSSP database: fold classification based on structurestructure alignment of proteins.
Nucleic Acids Research 1996, 24:206209. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Hobohm U, Scharf M, Schneider , Sander C: Selection of representative protein data sets.
Protein Science 1992, 1:409417. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Dietmann S, Park J, Notredame C, heger A, Lappe M, Holm L: A fully automatic evolutionary classification of protein folds: Dali Domain Dictionary version 3.
Nucleic Acids Research 2001, 29:5557. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Holm L, Sander C: Parser for protein folding units.
PROTEINS: Structure, Function and Genetics 1994, 19:256268. Publisher Full Text

Gewehr JE, Hintermair V, Zimmer R: AutoSCOP: automated prediction of SCOP classifications using unique patternclass mappings.
Bioinformatics 2007, 23:12031210. PubMed Abstract  Publisher Full Text

Kim YJ, Patel JM: A framework for protein structure classification and identification of novel protein structures.
BMC Bioinformatics 2006., 7(456) PubMed Abstract

Cheek S, Qi Y, Krishna SS, Kinch LN, Grishin NV: SCOPmap: Automated assignment of protein structures to evolutionary superfamilies.
BMC Bioinformatics 2004., 5 PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Rogen P, Fain B: Automatic classification of protein structure by using Gauss integral.
Proceedinf of national Academy of Science USA 2003, 100:119124. Publisher Full Text

Getz G, Vendruscolo M, Sachs D, Domany E: Automated Assignment of SCOP and CATH Protein Structure Classifications from FSSP.
Proteins: Structure, Function and Genetics 2002, 46:405415. Publisher Full Text

Gerstein M, Levitt M: Comprehensive assessment of automatic structural alignment against a manual standard, the SCOP classification of proteins.
Protein Science 1998, 7:445456. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Hadley C, Jones D: A systematic comparison of protein structure classifications: SCOP, CATH and FSSP.
Structures 1999, 7:10991112. Publisher Full Text

Harrison A, Pearl F, Mott R, Thornton J, Orengo C: Quantifying the similarities within fold space.
Journal of Molecular Biology 2002, 323:909926. PubMed Abstract  Publisher Full Text

Kolodny R, Koehl P, Levitt M: Comprehensive evaluation of protein structure alignment methods: scoring by geometric measures.
Journal of Molecular Biology 2005, 346:11731188. PubMed Abstract  Publisher Full Text

Sam V, Tai CH, Garnier J, Gibrat JF, Lee B, Munson PJ: ROC and confusion analysis of structure comparison methods identify the main causes of divergence from manual protein classification.
BMC Bioinformatics 2006. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Gibrat JF, Madej T, Bryant SS: Surprising similarities in structure comparison.
Current Opinion in Structural Biology 1996, 6:377385. PubMed Abstract  Publisher Full Text

Jung J, Lee B: Protein structure alignment using environmental profiles.
Protein Engineering 2000, 13:535543. PubMed Abstract  Publisher Full Text

Holm L, Sander C: Mapping the protein universe.
Science 1996, 273:595602. PubMed Abstract  Publisher Full Text

Holm L, Park J: DaliLite workbench for structure comparison.
Bioinformatics 2000, 16:566567. PubMed Abstract  Publisher Full Text

Frey B, Dueck D: Clustering by passing messages between data points.
Science 2007, 315:972976. PubMed Abstract  Publisher Full Text

Slonim N, Atwal GS, Gasper, Bialek W: Informationbased clustering.
PNAS 2005, 102:1829718302. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Nagano N, Orengo CA, Thornton JM: One fold with many functions: the evolutionary relationships between TIM barrel families based on their sequences, structures and functions.
Journal of Molecular Biology 2002, 321:741765. PubMed Abstract  Publisher Full Text

Orengo C, Michie AD, Jones S, Jones DT, Swindells MB, Thornton JM: CATHa hierarchic classification of protein domains structures.
Structures 1997, 5:10931108. Publisher Full Text

Lupas AN, Ponting CP, Russell RB: On the evolution of protein Folds: are similar motifs in different protein Folds the result of convergence, insertion, or Relics of an ancient peptide world?
Journal of Structural Biology 2001, 134:191203. PubMed Abstract  Publisher Full Text

MarchlerBauer A, Anderson JB, Cherukuri PF, DeWeeseScott C, Ly G, Gwadz M, He S, Hurwitz DI, Jackson JD, Ke Z, Lanczycki CJ, Liebert CA, Liu C, Lu F, Marchler GH, Mullokandov M, Shoemaker BA, Simonyan V, JS JSS, Thiessen PA, Yamashita RA, Yin JJ, Zhang D, Bryant SH: CDD: a Conserved Domain Database for protein classification.
Nucleic Acids Research 2005, 33:D1926. Publisher Full Text

Sonnhammer ELL, Kahn D: Modular arrangement of proteins as inferred from analysis of homology.
Protein Science 1994, 3:482492. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Sonnhammer ELL, Eddy SR, Durbin R: Pfam: a comprehensive database of protein domain families based on seed alignments.
PROTEINS: Structure, Function and Genetics 1997, 28:405420. Publisher Full Text

Chothia C, Lesk AM: The relation between the divergence of sequence and structure in proteins.
The EMBO journal 1986, 5:823826. PubMed Abstract  PubMed Central Full Text

Chandonia JM, Hon G, Walker NS, Conte LL, Koehl P, Levitt M, Brenner SE: The ASTRAL compendium in 2004.
Nucleic Acids Research 2004, 32:D189D192. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Biowulf cluster [http://biowulf.nih.gov/] webcite

MATLAB: Statistics toolbox, Linkage.
2006., R2006a

Shannon CE: A mathematical theory of communication.
The Bell System Technical Journal 1948, 27:379423, 623656, July, October.

MATLAB: Statistics toolbox, Classical Multidimensional Scaling.
2006., R2006a

Pettersen EF, Goddard TD, Huang CC, Couch GS, Greenblatt DM, Meng EC, Ferrin TE: UCSF Chimera  A Visualization System for Exploratory Research and Analysis.
Journal of Computational Chemistry 2004, 25:16051612. PubMed Abstract  Publisher Full Text