whatIsHMM_seanEddy_nbt04 - © 2004 Nature Publishing Group

Info iconThis preview shows page 1. Sign up to view the full content.

View Full Document Right Arrow Icon
This is the end of the preview. Sign up to access the rest of the document.

Unformatted text preview: © 2004 Nature Publishing Group http://www.nature.com/naturebiotechnology _computational BIOLOGY PRIMER What is a hidden Markov model? Sean R Eddy Statistical models called hidden Markov models are a recurring theme in computational biology. What are hidden Markov models, and why are they so useful for so many different problems? Often, biological sequence analysis is just a matter of putting the right label on each residue. In gene identification, we want to label nucleotides as exons, introns, or intergenic sequence. In sequence alignment, we want to associate residues in a query sequence with homologous residues in a target database sequence. We can always write an ad hoc program for any given problem, but the same frustrating issues will always recur. One is that we want to incorporate heterogeneous sources of information. A genefinder, for instance, ought to combine splice-site consensus, codon bias, exon/ intron length preferences and open reading frame analysis into one scoring system. How should these parameters be set? How should different kinds of information be weighted? A second issue is to interpret results probabilistically. Finding a best scoring answer is one thing, but what does the score mean, and how confident are we that the best scoring answer is correct? A third issue is extensibility. The moment we perfect our ad hoc genefinder, we wish we had also modeled translational initiation consensus, alternative splicing and a polyadenylation signal. Too often, piling more reality onto a fragile ad hoc program makes it collapse under its own weight. Hidden Markov models (HMMs) are a formal foundation for making probabilistic models of linear sequence ‘labeling’ problems1,2. They provide a conceptual toolkit for building complex models just by draw- Sean R. Eddy is at Howard Hughes Medical Institute & Department of Genetics, Washington University School of Medicine, 4444 Forest Park Blvd., Box 8510, Saint Louis, Missouri 63108, USA. e-mail: eddy@genetics.wustl.edu A = 0.25 C = 0.25 G = 0.25 T = 0.25 Start 1.0 E 0.1 A = 0.05 C=0 G = 0.95 T =0 5 A = 0.4 C = 0.1 G = 0.1 T = 0.4 1.0 0.9 Sequence: State path: End I 0.1 0.9 C T T C A TG TG A A AG C AG AC G T A A G T C A EE EEEEEEEEEEEEEEEE 5 I I I I I I I log P –41.22 –43.90 –43.45 –43.94 –42.58 –41.71 Parsing: 46% Posterior decoding: 28% 11% Figure 1 A toy HMM for 5′ splice site recognition. See text for explanation. ing an intuitive picture. They are at the heart of a diverse range of programs, including genefinding, profile searches, multiple sequence alignment and regulatory site identification. HMMs are the Legos of computational sequence analysis. A toy HMM: 5′ splice site recognition As a simple example, imagine the following caricature of a 5′ splice-site recognition problem. Assume we are given a DNA sequence that begins in an exon, contains one 5′ splice site and ends in an intron. The problem is to identify where the switch from exon to intron occurred—where the 5′ splice site (5′SS) is. For us to guess intelligently, the sequences of exons, splice sites and introns must have NATURE BIOTECHNOLOGY VOLUME 22 NUMBER 10 OCTOBER 2004 different statistical properties. Let’s imagine some simple differences: say that exons have a uniform base composition on average (25% each base), introns are A/T rich (say, 40% each for A/T, 10% each for C/G), and the 5′SS consensus nucleotide is almost always a G (say, 95% G and 5% A). Starting from this information, we can draw an HMM (Fig. 1). The HMM invokes three states, one for each of the three labels we might assign to a nucleotide: E (exon), 5 (5′SS) and I (intron). Each state has its own emission probabilities (shown above the states), which model the base composition of exons, introns and the consensus G at the 5′SS. Each state also has transition probabilities (arrows), the probabilities of moving from this state to a new state. The transition 1315 PRIMER © 2004 Nature Publishing Group http://www.nature.com/naturebiotechnology probabilities describe the linear order in which we expect the states to occur: one or more Es, one 5, one or more Is. So, what’s hidden? It’s useful to imagine an HMM generating a sequence. When we visit a state, we emit a residue from the state’s emission probability distribution. Then, we choose which state to visit next according to the state’s transition probability distribution. The model thus generates two strings of information. One is the underlying state path (the labels), as we transition from state to state. The other is the observed sequence (the DNA), each residue being emitted from one state in the state path. The state path is a Markov chain, meaning that what state we go to next depends only on what state we’re in. Since we’re only given the observed sequence, this underlying state path is hidden—these are the residue labels that we’d like to infer. The state path is a hidden Markov chain. The probability P (S,π|HMM,θ) that an HMM with parameters θ generates a state path π and an observed sequence S is the product of all the emission probabilities and transition probabilities that were used. For example, consider the 26-nucleotide sequence and state path in the middle of Figure 1, where there are 27 transitions and 26 emissions to tote up. Multiply all 53 probabilities together (and take the log, since these are small numbers) and you’ll calculate log P (S,π|HMM,θ) = –41.22. An HMM is a full probabilistic model — the model parameters and the overall sequence ‘scores’ are all probabilities. Therefore, we can use Bayesian probability theory to manipulate these numbers in standard, powerful ways, including optimizing parameters and interpreting the significance of scores. Finding the best state path In an analysis problem, we’re given a sequence, and we want to infer the hidden state path. There are potentially many state paths that could generate the same sequence. We want to find the one with the highest probability. For example, if we were given the HMM and the 26-nucleotide sequence in Figure 1, there are 14 possible paths that have nonzero probability, since the 5′SS must fall on one of 14 internal As or Gs. Figure 1 enumerates the six highest-scoring paths (those 1316 with G at the 5′SS). The best one has a log probability of –41.22, which infers that the most likely 5′SS position is at the fifth G. For most problems, there are so many possible state sequences that we could not afford to enumerate them. The efficient Viterbi algorithm is guaranteed to find the most probable state path given a sequence and an HMM. The Viterbi algorithm is a dynamic programming algorithm quite similar to those used for standard sequence alignment. Beyond best scoring alignments Figure 1 shows that one alternative state path differs only slightly in score from putting the 5′SS at the fifth G (log probabilities of –41.71 versus –41.22). How confident are we that the fifth G is the right choice? This is an example of an advantage of probabilistic modeling: we can calculate our confidence directly. The probability that residue i was emitted by state k is the sum of the probabilities of all the state paths that use state k to generate residue i (that is, πi = k in the state path π), normalized by the sum over all possible state paths. In our toy model, this is just one state path in the numerator and a sum over 14 state paths in the denominator. We get a probability of 46% that the best-scoring fifth G is correct and 28% that the sixth G position is correct (Fig. 1, bottom). This is called posterior decoding. For larger problems, posterior decoding uses two dynamic programming algorithms called Forward and Backward, which are essentially like Viterbi, but they sum over possible paths instead of choosing the best. Making more realistic models Making an HMM means specifying four things: (i) the symbol alphabet, K different symbols (e.g., ACGT, K = 4); (ii) the number of states in the model, M; (iii) emission probabilities ei(x) for each state i, that sum to one over K symbols x, Σ ei(x) = 1; and x (iv) transition probabilities ti( j ) for each state i going to any other state j (including itself) that sum to one over the M states j, Σ ti( j ) = 1. For example, in our toy splice-site model, maybe we’re not happy with our discrimination power; maybe we want to add a more realistic six-nucleotide consensus GTRAGT at the 5′ splice site. We can put a row of six HMM states in place of ‘5’ state, to model a six-base ungapped consensus motif, parameterizing the emission probabilities on known 5′ splice sites. And maybe we want to model a complete intron, including a 3′ splice site; we just add a row of states for the 3′SS consensus, and add a 3′ exon state to let the observed sequence end in an exon instead of an intron. Then maybe we want to build a complete gene model…whatever we add, it’s just a matter of drawing what we want. The catch HMMs don’t deal well with correlations between residues, because they assume that each residue depends only on one underlying state. An example where HMMs are usually inappropriate is RNA secondary structure analysis. Conserved RNA base pairs induce long-range pairwise correlations; one position might be any residue, but the base-paired partner must be complementary. An HMM state path has no way of ‘remembering’ what a distant state generated. Sometimes, one can bend the rules of HMMs without breaking the algorithms. For instance, in genefinding, one wants to emit a correlated triplet codon instead of three independent residues; HMM algorithms can readily be extended to tripletemitting states. However, the basic HMM toolkit can only be stretched so far. Beyond HMMs, there are more powerful (though less efficient) classes of probabilistic models for sequence analysis. 1. Rabiner, L.R. A tutorial on hidden Markov models and selected applications in speech recognition. Proc. IEEE 77, 257–286 (1989). 2. Durbin, R., Eddy, S.R., Krogh, A. & Mitchison, G.J. Biological Sequence Analysis: Probabilistic Models of Proteins and Nucleic Acids (Cambridge University Press, Cambridge UK, 1998). j Any model that has these properties is an HMM. This means that one can make a new HMM just by drawing a picture corresponding to the problem at hand, like Figure 1. This graphical simplicity lets one focus clearly on the biological definition of a problem. Wondering how some other mathematical technique really works? Send suggestions for future primers to askthegeek@natureny.com. VOLUME 22 NUMBER 10 OCTOBER 2004 N ATURE BIOTECHNOLOGY ...
View Full Document

This note was uploaded on 04/06/2010 for the course COMPUTER S COSC1520 taught by Professor Paul during the Spring '09 term at York University.

Ask a homework question - tutors are online