Title: Protein Clustering to Assemble Families of Homeomorphic Proteins
1Protein Clustering to Assemble Families of
Homeomorphic Proteins
- Chris Elsik
- c-elsik_at_tamu.edu
2Outline
- Using evolution to infer protein function
- The problem of automatic assembling of
homeomorphic (identical domain organization)
protein families - Agglomerative clustering
- Delineating protein domain boundaries
3Evolution Allows us to Infer Function
- The most powerful method for inferring function
of a gene or protein is by similarity searching a
sequence database. - Our ability to characterize biological properties
of a protein using sequence data alone stems from
properties conserved through evolutionary time. - Homologous (evolutionarily related) proteins
always share a common 3-dimensional folding
structure. - They often contain common active sites or binding
domains. - They frequently share common functions.
- Predictions made using similar, but
non-homologous proteins are much less reliable.
4Orthologs
- Homologs genes that are evolutionarily related
- There are two kinds of homologs
- Orthologs genes in different species that have
diverged from a common gene in an ancestral
species. - Paralogs genes that have diverged due to gene
duplication. - Orthologs are more likely than paralogs to have
conserved function. - Orthologs cannot be identified using BLAST or
FASTA sequence comparison alone. - Reliable ortholog identification requires
phylogenetic methods.
5Example Gene Tree (with plant genes)
Rice-2b
paralogs
Rice-2a
Maize-2
paralogs
Wheat-2
Sorghum-2
Barley-1
Wheat-1
orthologs
Maize-1
Sorghum-1
Arabidopsis
6Why shouldnt we depend on inferences based on
paralogs?
- Paralogs emerge after a gene duplication.
- Possible fates of duplicated genes
- Loss of function for one of the duplicates - lack
of selective pressure allows gene to mutate
beyond recognition - Emergence of new functional paralogs - one
duplicate aquires a new function, so selection
favors its maintenance in the genome - Sub-functionalization - both duplicates are
required to maintain the function of the original
7How do people identify orthologs?
- Best Hit
- Problems - 1) What if the ortholog is not present
in the database - there is always a best hit. 2)
For multidomain proteins, the best hit may be a
local domain hit. - Reciprocal Best Hit
- Problem - this is usually based on E-value, but
E-value is affected by the length of the match. - Predicted genes are often gene fragments (not the
true length of the gene) - This method is a phenetic, not phylogenetic,
approach. It does not take advantage of a model
of protein evolution. - Synteny - infer orthology by comparing locations
on chromosomes - Can be applied to closely related species, such
as mammals, but not distantly related species - difficult to distinguish tandemly duplicated
genes - Phylogenetics - build gene trees to distinguish
orthologs and paralogs
8 1gtgtgtCG2830-PA - 648 aa vs fly_test.lseg
libraryThe best scores are
length bits E(18717) _id
alenCG7235-PC ( 576) 474.6 9.4e-134
0.590 581CG7235-PA ( 576) 474.6
9.4e-134 0.590 581CG7235-PB ( 576)
474.6 9.4e-134 0.590 581 CG12101-PB
( 573) 474.4 1.1e-133 0.620
597CG12101-PA ( 573) 474.4 1.1e-133
0.620 597CG2830-PA-short ( 240) 339.2
2.2e-93 1.000 240CG16954-PB ( 558)
298.3 1.1e-80 0.450 536 CG16954-PA
( 558) 298.3 1.1e-80 0.450 536
Reciprocal Best Hit incorrect gene length
A gene fragment was simulated by replacing a
protein in the library dataset with a truncated
sequence. The library was searched with the
full-length protein. The best match is a paralog.
The identical, but truncated, match has a less
significant E-value. We would not be able to
identify the ortholog using the reciprocal best
hit method.
9Identifying Orthologs Using Phylogenetics
- Challenge - we need to create gene families prior
to creating phylogenetic gene trees - One approach is automated sequence clustering
10Reasons for Clustering Protein Sequences
- Classifying proteins - structural and functional
annotation - Identifying errors in gene prediction
- Selecting representative proteins
- Creating models of protein domain families to aid
identification of distantly related proteins - Studying protein family evolution
11Protein Clustering - General Method
- Perform a protein similarity search of a protein
database against itself using FASTA or BLAST - Use E-value or score as a similarity measure for
automated clustering - E-value is not the optimal linkage criterion due
to the length effect
12Problems in Clustering Proteins
- Multidomain proteins can cause unrelated proteins
to group together - This problem is amplified in large proteomes
- Cluster validation is difficult -
- We do not know the correct number of clusters
- We do not have a good test set for large proteomes
13Src Tyrosine Kinase
MGSNKSKPKDASQRRRSLEPAENVHGAGGGAFPASQTPSKPASADGHRGP
SAAFAPAAAEPKLFGGFNSS DTVTSPQRAGPLAGGVTTFVALYDYESRT
ETDLSFKKGERLQIVNNTEGDWWLAHSLSTGQTGYIPSNYV APSDSIQA
EEWYFGKITRRESERLLLNAENPRGTFLVRESETTKGAYCLSVSDFDNAK
GLNVKHYKIRKL DSGGFYITSRTQFNSLQQLVAYYSKHADGLCHRLTTV
CPTSKPQTQGLAKDAWEIPRESLRLEVKLGQGC FGEVWMGTWNGTTRVA
IKTLKPGTMSPEAFLQEAQVMKKLRHEKLVQLYAVVSEEPIYIVTEYMSK
GSLL DFLKGETGKYLRLPQLVDMAAQIASGMAYVERMNYVHRDLRAANI
LVGENLVCKVADFGLARLIEDNEYT ARQGAKFPIKWTAPEAALYGRFTI
KSDVWSFGILLTELTTKGRVPYPGMVNREVLDQVERGYRMPCPPEC PES
LHDLMCQCWRKEPEERPTFEYLQAFLEDYFTSTEPQYQPGENL
14Agglomerative Clustering Methods
- Single linkage - similarity between closest
neighbors meets a threshold (BLASTCLUST) - Complete linkage - similarity between furthest
neighbors meets a threshold - Average linkage - average pairwise similarity
meets a threshold - Fractional linkage - fraction of pairwise
similarities meet a threshold
15Single linkage
If A is related to B, and B is related to C,
then A is related to C.
This is true only if proteins are homologous over
their entire length. Multidomain proteins pull
unrelated proteins into the same cluster
A
B
C
16Complete linkage
Each protein must be homologous to every other
protein in the cluster.
In the analysis of the human genome, Celera used
the criterion that the list of matches for each
protein must be equivalent.
Query A B
C D Matches
C D
A B C
A B C D
A B C
Clusters
A B
C D
Problem Does the domain structure of D differ
from A and B, or is it homologous, but too
distantly related to be detected by the BLAST
search? Complete linkage generates many small
clusters.
17Fractional linkage
Proteins must have a fraction of matches in
common.
In analysis of the human genome, Celera used the
Lek metric with a cutoff of .75 Lek 2 x
(MatchesA ? MatchesB) / (MatchesA MatchesB)
LekAB 2 x 3 / (3 3) 1
LekAC 2 x 3 / (3 4) 0.85
Clusters
D
A B C
LekCD 2 x 2 / (4 2) 0.67
LekAD 2 x 1 / (32) 0.40
18Average linkage
The average pairwise similarity between all
proteins must meet a threshold. In hierarchical
average linkage, the most closely related
proteins are grouped first. Then the threshold is
relaxed, allowing clusters to merge.
19Which is method is best for clustering
full-length protein sequences?
- Objectives
- Group full-length proteins with similar domain
organization - Minimize the problem caused by multidomain
proteins - Minimize number of singletons
20Comparing Protein Clustering Methods
- Four methods - single, complete, average and
fractional linkage - Thresholds range from E ? 10-10 (permissive) to E
? 10-200 (stringent) - Look at cluster number, number of singletons,
largest cluster size, CDC score for each cluster
set - Test set the Arabidopsis proteome (25,000
proteins)
21Comparison of cluster sets assembled using
linkage threshold E ? 10-10
Singletons
Largest Cluster Size
Method
Clusters (including Singletons)
SL
7702
4940
5852
9607
AL
385
5943
FL
12072
8248
428
60
CL
17639
14151
22Comparing protein cluster sets CDC Score
- Cluster Domain Consistency (CDC) score - a metric
for testing similarity in domain organization
among clustered proteins - CDC 0 proteins have no domains in common
- CDC 1 proteins have identical domain
organization - The proteins first must be searched against a
domain database (Pfam) to identify domains
23(No Transcript)
24Pfam CDC Score
25Linkage Criteria for Protein Clustering
- The linkage criterion should discern full- and
partial-length matches - If proteins are similar in domain organization,
we expect them to match over their entire length
26Linkage Criteria E-value
- E-value may not be the best linkage criterion
because it is length-dependent - Matches to longer proteins have lower E-values
(more significant) - Length dependency makes it difficult to discern
partial and full-length protein matches
27Linkage Criteria Global Alignment Score
- Partial matches will have low scores and
full-length matches will have high scores - Disadvantage global alignments for each protein
pair must be generated using the ALIGN program
28Linkage Criteria Weighted Local Alignment Score
- The local Smith-Waterman score for a pair of
similar proteins is divided by the self-match
score for the larger of the two proteins. - The weighted score is low for partial matches.
- The advantage over global score is that global
alignments are not required.
29Linkage Criteria Comparison
- Two clustering methods single and average
linkage - Three linkage criteria E()-value, weighted local
score, global score - Range of thresholds tested for each criterion
- Look at cluster number, number of singletons,
largest cluster size, CDC score for each cluster
set
30Cluster set comparison using most permissive
thresholds that resulted in biologically
meaningful family sizes (lt 600 members).
Method
Single Linkage
Average Linkage
Weighted Local Score
Weighted Local Score
Global Score
Linkage Parameter
Global Score
E()-value
E()-value
E() ? 10-10
E() ? 10-40
0.05
0.5
Threshold
1.1
0.25
Clusters (including singletons)
14570
9607
11469
8087
10073
10905
4643
11560
6195
5943
Singletons
7780
7429
250
538
Largest Cluster
385
225
171
563
31Pfam CDC Score
32Protein Clustering - Conclusions
- Single linkage using E-value is not appropriate
for clustering eukaryote proteomes. - Weighted local score and global score are better
linkage criteria than E-value. - Single linkage using a weighted local score above
0.40 performs as well as average linkage and
generates a moderate number (15,746) of
Arabidopsis clusters. - Average linkage generates the smallest cluster
numbers and reasonable cluster sizes, and should
be the useful for grouping distantly related
proteins for function and structure prediction.
33Which clustering method is best prior to
phylogenetic analysis of predicted proteins?
- Complete linkage? - Too many singletons
- Fractional linkage? - Also a large number of
singletons - Average linkage?
- Problem caused by truncated or overextended gene
predictions. Average linkage would likely
separate orthologs with incorrectly predicted
protein lengths. - Single linkage with a strict alignment criterion?
- Also has problem caused by truncated or
overextended gene predictions - A strict alignment criterion separates orthologs
due to incorrectly predicted protein lengths - Single linkage without a strict alignment
criterion - Problem associated with multidomain proteins
34SL-AL a compromise
- We have developed the SL-AL algorithm, which uses
a combination of single linkage and average
linkage - First sequences are clustered using single
linkage and a permissive alignment criterion - Percent identity for each protein pair in a
single linkage cluster is recorded - Any cluster with lt a threshold percent identity
is reclustered using average linkage
35Advantages of SL-AL
- The multidomain problem is reduced.
- The algorithm is more tolerant to incorrect gene
predictions, and less likely to separate
incorrectly predicted orthologs into separate
clusters. The cluster set can be used to identify
gene prediction problems. - A minimum pairwise percent identity can be set,
making multiple alignment possible.
36Another Alternative for Protein Clustering
- Overcome the multidomain problem by dissecting
proteins into domains prior to clustering - Prodom uses protein alignments to identify domain
families (129,000 domain families) - Problem - partial sequences in protein databases
cause the overfragmenting of domains. - Our work toward a solution use an alternative
approach to identify domain boundaries
37Identifying protein linkers to delineate domain
boundaries
- Protein structural domain - a unit of a protein
that can independently fold into a stable
tertiary structure - Evolutionary domain - a conserved modular unit
that can be identified by aligning protein
sequences - Often structural domains are also evolutionary
domains (not always) - Linker - a peptide that connects two protein
domains
38Src Tyrosine Kinase
MGSNKSKPKDASQRRRSLEPAENVHGAGGGAFPASQTPSKPASADGHRGP
SAAFAPAAAEPKLFGGFNSS DTVTSPQRAGPLAGGVTTFVALYDYESRT
ETDLSFKKGERLQIVNNTEGDWWLAHSLSTGQTGYIPSNYV APSDSIQA
EEWYFGKITRRESERLLLNAENPRGTFLVRESETTKGAYCLSVSDFDNAK
GLNVKHYKIRKL DSGGFYITSRTQFNSLQQLVAYYSKHADGLCHRLTTV
CPTSKPQTQGLAKDAWEIPRESLRLEVKLGQGC FGEVWMGTWNGTTRVA
IKTLKPGTMSPEAFLQEAQVMKKLRHEKLVQLYAVVSEEPIYIVTEYMSK
GSLL DFLKGETGKYLRLPQLVDMAAQIASGMAYVERMNYVHRDLRAANI
LVGENLVCKVADFGLARLIEDNEYT ARQGAKFPIKWTAPEAALYGRFTI
KSDVWSFGILLTELTTKGRVPYPGMVNREVLDQVERGYRMPCPPEC PES
LHDLMCQCWRKEPEERPTFEYLQAFLEDYFTSTEPQYQPGENL
39Amino acid propensity
- Our goal is to predict protein linker regions
using amino acid sequence alone, in the absence
of known homologs - We take advantage of amino acid propensity - the
preference of some amino acids for linkers or
domains
40Our Dataset
- Pfam-A - a collection of domain families created
using profile HMMs built on multiple alignments
of homologous proteins. (Boeckmann et al. Nucleic
Acids Research 2003) - Linker a sequence segment of 4 to 20 residues
that connects two adjacent regions identified by
Pfam as domains. - Non-linker regions sequence segments excluding
linker regions. - Used only protein sequences whose entire length
can be classified as linker or domain by our
criteria, except we allowed up to 20 non-domain
residues at the N- and C-termini. - Results - 11968 sequences with at least one
linker region (14339 linker, 28726 corresponding
domains and 824 unique domain regions).
41Removing Redundancy in the Dataset
- 11968 proteins were first grouped into
homeomorphic families (identical domain
organization). - All by all sequence comparison of the 11968
sequences using FASTA - Single-linkage clustering using criteria of
E-value10-6 and at least 80 alignment coverage. - Some of the resulting clusters contained
sequences with different domain organizations,
due to the transitive nature of single-linkage
clustering, so we selected one sequence from each
domain organization within each cluster. - Result - 802 sequences with at least one linker
region. These 802 sequences contained 993 linkers
and 1988 corresponding domain regions from 376
unique Pfam-A domain families. - The average length of linkers and domains was
11.24 and 141.38, respectively.
42Linker Index
- Linker index, yl, reflects the preference of
amino acids in the linkers relative to the domain
region from Suyama and Ohara, 2003,
Bioinformatics. - yl - ln ( fllinker / fldomain )
- Where fllinker relative frequency of the amino
acid l in the linker (region in the data set. - yl will be negative if the relative frequency of
amino acid l in the linker region is greater than
its relative frequency in the domain region. - To calculate the smoothed linker index, average
the linker index within each window size w and
assign this averaged linker index value y to the
center amino acid of the window by sliding from
the N-terminus to the C-terminus of a protein
sequence. - Window size, w 9, provided the greatest
discrimination between linker and non-linker
regions among the window sizes from 3 to 20.
43Relative Amino Acid Frequency
Amino Linker Domain yl Acid
A 7.97 (7.94) 8.10 0.0166 C 0.89
(1.24) 1.50 0.5724 D 6.32
(5.28) 5.60 -0.1278 E 7.97
(6.89) 6.60 -0.1794 F 2.74
(4.34) 4.03 0.3561 G 7.74
(6.14) 7.37 -0.0442 H 1.91
(2.32) 2.27 0.1643 I 4.73
(5.13) 6.37 0.2758 K 6.97
(5.72) 5.81 -0.2134 L 7.51
(9.60) 9.54 0.2523 M 2.13
(2.15) 2.24 0.0197 N 4.22
(4.12) 4.08 -0.0786 P 6.63
(6.07) 4.30 -0.4188 Q 3.90
(4.05) 3.33 -0.1051 R 5.77
(5.79) 5.24 -0.0762 S 7.20
(5.55) 6.13 -0.1629 T 5.80
(5.66) 5.35 -0.0701 V 6.24
(6.64) 7.34 0.1782 W 0.81
(1.24) 1.32 0.3836 Y 2.46
(3.47) 3.38 0.2500
Numbers in parentheses are from the linker
database of George and Heringa, 2002, Protein
Engineering.
44Hidden Markov Model to predict linker residues
- Sequences are assumed to have a structure
composed of regions that are homogeneous within a
region but may differ between regions. - We assume protein sequence data is produced by a
hidden Markov model and compositional variation
is likely to reflect functional or structural
differences between regions. - Each region is classified into one of a finite
number of states (linker and non-linker) we wish
to estimate the states given the observed protein
sequence. - Instead of recognizing the protein sequence as a
string of amino acids (categorical variables), we
recognize the protein sequence as a string of
linker index values (continuous variables).
45Predicting linker state
- Our objective is to identify linker index values
that discriminate linker and non-linker regions. - Used Gibbs sampling to overcome the problem of
missing data. - Calculate the probability state (linker or
non-linker) for each residue i in a protein
sequence. - Predict the state of an amino acid using the
classification variable - CV 1 if the probability of state k is ? x.
- CV 0 if the probability of state k is ? x.
46Training and Testing
- We applied our model to the protein sequence
dataset constructed from Pfam-A. - 5-fold cross validation was applied to the
dataset - we divided the dataset into the
training dataset and the test dataset randomly
with the ratio of 41. - We trained the model with the training dataset of
642 sequences and tested the trained model with
the test dataset of 160 sequences. - We ran Gibbs sampling with 40,000 iterations and
10,000 burn-in to train the model.
47Examples of good predictions
48Examples of Over-prediction
49Sensitivity and Selectivity
- Sensitivity the percentage of actual linker
residues that were predicted to be linker.
SnTp(TpFn) - Specificity the percentage of predicted linker
residues that were truly linker. SpTp/(TpFp)
50Sensitivity and Selectivity plotted against
cut-off for classification variable
Sensitivity and selectivity are each 67 at their
intersection.
51Conclusions
- This method will be useful in defining linkers
for evolutionary domains. - Our method appears to do well at distinguishing
linkers from intra-domain loops. - We must test our approach using a structurally
defined dataset to fully understand its ability
to distinguish structural linkers from
intra-domain loops. - Since Pfam-A identifies domains as evolutionarily
conserved units, non-conserved intra-domain loops
can cause structural domains to be annotated as
multiple Pfam-A domains. Thus, some of our Pfam-A
defined linkers may actually be loops in
structural domains. - Conversely, two structural domains that are
always found together may be defined by Pfam-A as
a single evolutionary domain some of our false
positives may actually be structural linkers.
52Acknowledgements
- William Pearson - University of Virginia
- Bani Mallick - Texas AM University
- Elsik Lab
- Kyounghwa Bae
- Justin Reese
- Anand Venkatraman
- Shreyas Murthi
- Michael Dickens
- Juan Anzola