Computational biology: Sequence alignment and profile HMMs10-601 Machine Learning2Central dogmaProteinmRNADNAtranscriptiontranslationCCTGAGCCAACTATTGATGAAPEPTIDECCUGAGCCAACUAUUGAUGAA3Central dogmaProteinmRNADNAtranscriptiontranslationCCTGAGCCAACTATTGATGAAPEPTIDECCUGAGCCAACUAUUGAUGAACan be measured using sequencing techniquesCan be measured using microarraysCan be measured using mass spectrometry4Comparison of Different Organisms Genome size Num. of genesE. coli .05*1084,200Yeast .15*1086,000Worm 1*10818,400Fly 1.8*10813,600Human 30*10825,000Plant 1.3*10825,0005Growth in biological dataLu et al Bioinformatics 20096Assigning function to proteins• One of the main goals of molecular (and computational) biology.• There are 25000 human genes and the vast majority of their functions is still unknown• Several ways to determine function- Direct experiments (knockout, overexpression)- Interacting partners- 3D structures- Sequence homologyHardEasier7Function from sequence homology• We have a query gene: ACTGGTGTACCGAT• Given a database containing genes with known function, our goal is to find similar genes from this database (possibly in another organism)• When we find such gene we predict the function of the query gene to be similar to the resulting database gene• Problems- How do we determine similarity?8Sequence analysis techniques• A major area of research within computational biology.• Initially, based on deterministic or heuristic alignment methods• More recently, based on probabilistic inference methods9Sequence analysis• Traditional- Dynamic programming• Probabilsitic- Profile HMMs10Pairwise sequence alignmentACATTGAACATTA C A T T GA A C A T TAGCCTTAGCATTA G C C T TA G C A T T11Pairwise sequence alignmentAGCCTTACCATTA G C C T TA C C A T TAGCCTTAGCATTA G C C T TA G C A T T• We cannot expect the alignments to be perfect.• Major reasons include insertion, deletion and substitutions.• We need to allow gaps in the resulting alignment.12Scoring AlignmentsjxixjiqqIyxP )|,(iyxiipMyxP )|,()log(),(,babaqqpbas • Alignments can be scored by comparing the resulting alignment to a background (random) model. Independent RelatedScore for alignment:),(iiiyxsSwhere:Can be computed for each pair of letters13Scoring AlignmentsjxixjiqqIyxP )|,(iyxiipMyxP )|,()log(),(,babaqqpbas • Alignments can be scored by comparing the resulting alignment to a background (random) model. Independent RelatedScore for alignment:),(iiiyxsSwhere:In other words, we are trying to find an alignment that maximizes the likelihood ratio of the aligned pair compared to the background model14Computing optimal alignment: The Needham-Wuncsh algorithmF(i-1,j-1) F(i-1,j)F(i,j-1) F(i,j)F(i,j) = maxF(i-1,j-1)+s(xi,xj)F(i-1,j)+dF(i,j-1)+dA G C C T TACCATTd is a penalty for a gap15ExampleA G C C T T0 -1 -2 -3 -4 -5 -6A -1C -2C -3A -4T -5T -6Assume a simple model where S(a,b) = 1 if a=b and -5 otherwise.Also, assume that d = -116ExampleAssume a simple model where S(a,b) = 1 if a=b and -5 otherwise.Also, assume that d = -1F(i,j) = maxF(i-1,j-1)+s(xi,xj)F(i-1,j)+dF(i,j-1)+dA G C C T T0 -1 -2 -3 -4 -5 -6A -1 1C -2C -3A -4T -5T -617ExampleAssume a simple model where S(a,b) = 1 if a=b and -5 otherwise.Also, assume that d = -1A G C C T T0 -1 -2 -3 -4 -5 -6A -1 1 0C -2 0C -3A -4T -5T -6F(i,j) = maxF(i-1,j-1)+s(xi,xj)F(i-1,j)+dF(i,j-1)+d18ExampleAssume a simple model where S(a,b) = 1 if a=b and -5 otherwise.Also, assume that d = -1A G C C T T0 -1 -2 -3 -4 -5 -6A -1 1 0 -1 -2 -3 -4C -2 0 -1C -3 -1A -4 -2T -5 -3T -6 -4F(i,j) = maxF(i-1,j-1)+s(xi,xj)F(i-1,j)+dF(i,j-1)+d19ExampleAssume a simple model where S(a,b) = 1 if a=b and -5 otherwise.Also, assume that d = -1A G C C T T0 -1 -2 -3 -4 -5 -6A -1 1 0 -1 -2 -3 -4C -2 0 -1 1 0 -1 -2C -3 -1 -2 0 2 1 0A -4 -2 -3 -1 1 0 -1T -5 -3 -4 -2 0 2 1T -6 -4 -5 -3 -1 1 320ExampleAssume a simple model where S(a,b) = 1 if a=b and -5 otherwise.Also, assume that d = -1A G C C T T0 -1 -2 -3 -4 -5 -6A -1 1 0 -1 -2 -3 -4C -2 0 -1 1 0 -1 -2C -3 -1 -2 0 2 1 0A -4 -2 -3 -1 1 0 -1T -5 -3 -4 -2 0 2 1T -6 -4 -5 -3 -1 1 321ExampleAssume a simple model where S(a,b) = 1 if a=b and -5 otherwise.Also, assume that d = -1A G C C T T0 -1 -2 -3 -4 -5 -6A -1 1 0 -1 -2 -3 -4C -2 0 -1 1 0 -1 -2C -3 -1 -2 0 2 1 0A -4 -2 -3 -1 1 0 -1T -5 -3 -4 -2 0 2 1T -6 -4 -5 -3 -1 1 3A G C C T TA C C A T T22Running time• The running time of an alignment algorithms if O(n2)• This doesn’t sound too bad, or is it?• The time requirement for doing global sequence alignment is too high in many cases.• Consider a database with tens of thousands of sequences. Looking through all these sequences for the best alignment is too time consuming.• In many cases, a much faster heuristic approach can achieve equally good results.23Sequence analysis• Traditional- Dynamic programming• Probabilsitic- Profile HMMs24Protein families• Proteins can be classified into families (and further into sub families etc.)• A specific family includes proteins with similar high level functions• For example:- Transcription factors- Receptors- Etc.Family assignment is an important first step towards function prediction25Methods for Characterizing a Protein Family• Objective: Given a number of related sequences, encapsulate what they have in common in such a way that we can recognize other members of the family. • Some standard methods for characterization:– Multiple Alignments– Regular Expressions– Consensus Sequences– Hidden Markov Models26Multiple Alignment Process• Process of aligning three or more sequences with each other• We can determine such alignment by generalizing the algorithm to align two sequences• Running time exponential in the number of sequences27Training a HMM from an existing alignment– Start with a predetermined number of states accounting for matches, insertions and deletions.– For each position in the model, assign a column in the multiple alignment that is relatively conserved.– Emission probabilities are set according to amino acid counts in columns.– Transition probabilities are set according to how many sequences make use of a given delete or insert state.MLE estimates28Remember the simple example• Chose six positions in model.• Highlighted area was selected to be modeled by an insert due to variability.• Can also do neat tricks for picking length of model,
View Full Document