WEBVTT 00:00:00.000 --> 00:00:00.610 align:middle line:90% 00:00:00.610 --> 00:00:02.440 align:middle line:84% The following content is provided 00:00:02.440 --> 00:00:06.340 align:middle line:84% by MIT OpenCourseWare under a Creative Commons license. 00:00:06.340 --> 00:00:10.120 align:middle line:84% Additional information about our license and MIT OpenCourseWare, 00:00:10.120 --> 00:00:12.730 align:middle line:84% in general, is available at ocw.mit.edu. 00:00:12.730 --> 00:00:17.270 align:middle line:90% 00:00:17.270 --> 00:00:26.130 align:middle line:84% GEORGE CHURCH: OK, so welcome to the fourth lecture. 00:00:26.130 --> 00:00:29.300 align:middle line:84% This will be the second one on the subject of DNA. 00:00:29.300 --> 00:00:34.370 align:middle line:84% The major difference between this lecture and the last one 00:00:34.370 --> 00:00:38.060 align:middle line:84% is that last one, so-called DNA 1, 00:00:38.060 --> 00:00:47.030 align:middle line:84% we focused in on types of mutants, their really closely 00:00:47.030 --> 00:00:48.480 align:middle line:90% related DNA sequences. 00:00:48.480 --> 00:00:50.960 align:middle line:84% And this one, we'll talk about the most distant related 00:00:50.960 --> 00:00:52.700 align:middle line:90% biopolymer sequences. 00:00:52.700 --> 00:00:56.810 align:middle line:84% So we went through the way you can generate 00:00:56.810 --> 00:01:01.550 align:middle line:84% closely related sequences and the way that the populations 00:01:01.550 --> 00:01:04.340 align:middle line:84% which are made up of closely related sequences 00:01:04.340 --> 00:01:09.470 align:middle line:84% obtain those allele frequencies through mutation drift 00:01:09.470 --> 00:01:10.340 align:middle line:90% and selection. 00:01:10.340 --> 00:01:14.900 align:middle line:84% And I argued that, deep down, at the most precise level, 00:01:14.900 --> 00:01:19.670 align:middle line:84% you've got a binomial distribution behind each 00:01:19.670 --> 00:01:24.562 align:middle line:84% of the processes of mutation drift and selection, at least 00:01:24.562 --> 00:01:26.020 align:middle line:84% under a certain set of assumptions. 00:01:26.020 --> 00:01:28.500 align:middle line:84% And then we went on to another very valuable statistic, 00:01:28.500 --> 00:01:30.750 align:middle line:84% which is the chi-square, which you can use in a number 00:01:30.750 --> 00:01:32.510 align:middle line:84% of different scenarios but, in particular, 00:01:32.510 --> 00:01:36.050 align:middle line:84% with association studies, where we did a simple case 00:01:36.050 --> 00:01:41.030 align:middle line:84% of a two-allele system with two outcomes-- 00:01:41.030 --> 00:01:44.150 align:middle line:90% HIV resistance and sensitivity. 00:01:44.150 --> 00:01:50.420 align:middle line:84% And then this association led to a broader discussion of alleles 00:01:50.420 --> 00:01:52.520 align:middle line:84% and haplotypes and genotypes, in general, 00:01:52.520 --> 00:01:55.640 align:middle line:84% and how one can obtain those and, in the process of doing 00:01:55.640 --> 00:01:58.700 align:middle line:84% so, expose oneself to random systematic errors, 00:01:58.700 --> 00:02:02.310 align:middle line:84% a theme that we'll return to from time to time. 00:02:02.310 --> 00:02:07.800 align:middle line:84% So today, we'll talk not about very closely related sequences 00:02:07.800 --> 00:02:09.289 align:middle line:84% but very distant related sequences 00:02:09.289 --> 00:02:10.789 align:middle line:84% and what are the algorithms that are 00:02:10.789 --> 00:02:12.470 align:middle line:90% available for finding those. 00:02:12.470 --> 00:02:17.420 align:middle line:84% These distant related sequences get a very different set 00:02:17.420 --> 00:02:18.170 align:middle line:90% of avenues. 00:02:18.170 --> 00:02:22.080 align:middle line:84% Here, we're trying to look for hints, for hypotheses in, say, 00:02:22.080 --> 00:02:27.860 align:middle line:84% new genome sequences as to what new genes might do. 00:02:27.860 --> 00:02:33.980 align:middle line:84% So we'll begin by comparing different types of algorithms 00:02:33.980 --> 00:02:36.860 align:middle line:84% that will allow us to do alignments. 00:02:36.860 --> 00:02:39.710 align:middle line:84% In particular, we'll stress the hero of today's show 00:02:39.710 --> 00:02:41.270 align:middle line:90% is dynamic programming. 00:02:41.270 --> 00:02:43.880 align:middle line:84% It comes up again and again throughout the course, 00:02:43.880 --> 00:02:47.150 align:middle line:84% not merely in pairwise sequence alignment 00:02:47.150 --> 00:02:49.970 align:middle line:90% but multisequence alignment. 00:02:49.970 --> 00:02:57.770 align:middle line:84% And we'll then go on to how the space in your computer, 00:02:57.770 --> 00:03:00.710 align:middle line:84% or the memory, dedicated to processing 00:03:00.710 --> 00:03:05.600 align:middle line:90% the time determine trade-offs. 00:03:05.600 --> 00:03:07.850 align:middle line:84% And in some cases, you will have to even sacrifice 00:03:07.850 --> 00:03:10.010 align:middle line:90% accuracy or completeness. 00:03:10.010 --> 00:03:13.940 align:middle line:84% And this also leads to the issue of finding 00:03:13.940 --> 00:03:18.950 align:middle line:84% genes, a particular type of distant sequence comparison. 00:03:18.950 --> 00:03:21.620 align:middle line:84% The ones interested in is finding motifs that 00:03:21.620 --> 00:03:23.900 align:middle line:90% are involved in finding genes. 00:03:23.900 --> 00:03:27.050 align:middle line:84% And finally, we'll end on hidden Markov model, the simplest one 00:03:27.050 --> 00:03:29.330 align:middle line:84% that I could think of that would really illustrate 00:03:29.330 --> 00:03:32.810 align:middle line:84% the idea of Markov models, probabilistic models, 00:03:32.810 --> 00:03:34.370 align:middle line:90% and the hiddenness. 00:03:34.370 --> 00:03:39.980 align:middle line:84% And so this is illustrated with a single dinucleotide 00:03:39.980 --> 00:03:42.590 align:middle line:90% with two states. 00:03:42.590 --> 00:03:44.900 align:middle line:90% So this puts in context-- 00:03:44.900 --> 00:03:48.710 align:middle line:84% in the first couple of lectures, we 00:03:48.710 --> 00:03:52.310 align:middle line:84% talked about the tree of life and how, right at the core, 00:03:52.310 --> 00:03:55.280 align:middle line:84% we had the common and most simple 00:03:55.280 --> 00:03:58.910 align:middle line:84% forms, which share the simple genetic code. 00:03:58.910 --> 00:04:01.640 align:middle line:84% And that was-- and then last time, 00:04:01.640 --> 00:04:05.810 align:middle line:84% we talked about the very tip of one branch of one 00:04:05.810 --> 00:04:08.360 align:middle line:84% of these trees-- basically, the human branch 00:04:08.360 --> 00:04:10.820 align:middle line:90% of the animal branch. 00:04:10.820 --> 00:04:15.590 align:middle line:84% And what has happened at the tip of that the last 5,000 00:04:15.590 --> 00:04:20.029 align:middle line:84% generations, since the fairly significant bottlenecks 00:04:20.029 --> 00:04:25.890 align:middle line:84% that resulted that predated the population explosion. 00:04:25.890 --> 00:04:29.420 align:middle line:84% And so this time, we're going to talk about the whole tree 00:04:29.420 --> 00:04:35.480 align:middle line:84% and how these very deep branches can be obtained in comparing 00:04:35.480 --> 00:04:37.940 align:middle line:90% biopolymer sequences. 00:04:37.940 --> 00:04:39.470 align:middle line:84% Some of the earliest trees of life 00:04:39.470 --> 00:04:41.295 align:middle line:84% were based on the ribosomal RNA sequence 00:04:41.295 --> 00:04:43.670 align:middle line:84% because that's something that's fairly easy to trace back 00:04:43.670 --> 00:04:45.050 align:middle line:90% to common ancestors. 00:04:45.050 --> 00:04:46.940 align:middle line:84% But you also want to be able to do this 00:04:46.940 --> 00:04:49.310 align:middle line:84% with a variety of proteins, some of which 00:04:49.310 --> 00:04:52.613 align:middle line:84% go only a short distance traceably in sequence, 00:04:52.613 --> 00:04:54.155 align:middle line:84% and some of them go all the way back. 00:04:54.155 --> 00:04:57.243 align:middle line:90% 00:04:57.243 --> 00:04:58.910 align:middle line:84% What can we do with dynamic programming? 00:04:58.910 --> 00:05:02.660 align:middle line:84% I alluded to its multiple uses throughout the course. 00:05:02.660 --> 00:05:04.560 align:middle line:90% Here's some examples. 00:05:04.560 --> 00:05:06.530 align:middle line:84% We already mentioned shotgun sequence assembly. 00:05:06.530 --> 00:05:08.990 align:middle line:84% This is relatively easy because the sequences are fairly 00:05:08.990 --> 00:05:11.660 align:middle line:90% closely related to one another. 00:05:11.660 --> 00:05:15.860 align:middle line:84% Multiple alignments include sequence assembly, 00:05:15.860 --> 00:05:17.617 align:middle line:84% where you have multiple similar sequences. 00:05:17.617 --> 00:05:19.700 align:middle line:84% But as you get to more and more distant sequences, 00:05:19.700 --> 00:05:24.110 align:middle line:84% you can glean more and more about structural and functional 00:05:24.110 --> 00:05:27.530 align:middle line:84% conclusions about proteins and nucleic acids 00:05:27.530 --> 00:05:30.350 align:middle line:84% when you have a very large number, in the hundreds. 00:05:30.350 --> 00:05:34.010 align:middle line:84% And we'll see the challenges that come there. 00:05:34.010 --> 00:05:37.220 align:middle line:84% Repeats are a particular-- you can 00:05:37.220 --> 00:05:39.530 align:middle line:84% have repeats within a genome, or we 00:05:39.530 --> 00:05:44.232 align:middle line:84% can have alignments between different species. 00:05:44.232 --> 00:05:46.190 align:middle line:84% Now, birdsong seems a little out of place here. 00:05:46.190 --> 00:05:49.065 align:middle line:84% But historically, it's one of the first applications 00:05:49.065 --> 00:05:50.940 align:middle line:84% of dynamic programming where you, in a sense, 00:05:50.940 --> 00:05:56.270 align:middle line:84% have a continuous time axis, and you sample that 00:05:56.270 --> 00:05:58.250 align:middle line:84% at discrete points, often, in order 00:05:58.250 --> 00:06:06.860 align:middle line:84% to sample the intensity of the audio recording. 00:06:06.860 --> 00:06:11.540 align:middle line:84% And this will have a more direct analog 00:06:11.540 --> 00:06:14.300 align:middle line:84% than sequence alignments, which we'll do today. 00:06:14.300 --> 00:06:17.810 align:middle line:84% Later, when we're doing the RNA analyses, 00:06:17.810 --> 00:06:21.980 align:middle line:84% gene expression will show a direct time warping 00:06:21.980 --> 00:06:26.615 align:middle line:84% algorithm, which is very similar in outline to the birdsong. 00:06:26.615 --> 00:06:29.300 align:middle line:90% 00:06:29.300 --> 00:06:32.030 align:middle line:84% And then, finally, hidden Markov models, 00:06:32.030 --> 00:06:33.950 align:middle line:84% which would be the last thing today, 00:06:33.950 --> 00:06:38.840 align:middle line:84% use dynamic programming as part of the process by which 00:06:38.840 --> 00:06:44.010 align:middle line:84% the decisions are made about the model that's represented, 00:06:44.010 --> 00:06:45.320 align:middle line:90% the hidden Markov model. 00:06:45.320 --> 00:06:48.700 align:middle line:84% And these hidden Markov models, in turn, 00:06:48.700 --> 00:06:51.890 align:middle line:84% allow us to do RNA gene searches, structure prediction, 00:06:51.890 --> 00:06:57.720 align:middle line:84% distant protein homologies, speech recognition, and so on. 00:06:57.720 --> 00:07:02.750 align:middle line:84% So there's three types of alignments and scores 00:07:02.750 --> 00:07:06.800 align:middle line:84% that we'll discuss on the slide number 5. 00:07:06.800 --> 00:07:09.290 align:middle line:84% The main dichotomy is between global and local. 00:07:09.290 --> 00:07:11.000 align:middle line:84% Originally, in sequence alignment 00:07:11.000 --> 00:07:14.338 align:middle line:84% in the '70s, Needleman-Wunsch had a global algorithm, 00:07:14.338 --> 00:07:15.380 align:middle line:90% which has been modified-- 00:07:15.380 --> 00:07:18.230 align:middle line:84% Smith and Waterman's local algorithm-- about a decade 00:07:18.230 --> 00:07:19.010 align:middle line:90% later. 00:07:19.010 --> 00:07:25.940 align:middle line:84% And the major difference here is that in a global algorithm, 00:07:25.940 --> 00:07:29.210 align:middle line:84% you have reason to believe that the sequences align end to end. 00:07:29.210 --> 00:07:32.210 align:middle line:84% And you really just want to ask how many mismatches 00:07:32.210 --> 00:07:35.750 align:middle line:84% and insertions and deletions are there 00:07:35.750 --> 00:07:37.320 align:middle line:90% in the middle of the sequence? 00:07:37.320 --> 00:07:40.430 align:middle line:84% So for example, you might have two proteins that really have 00:07:40.430 --> 00:07:42.140 align:middle line:90% the same start and stop site. 00:07:42.140 --> 00:07:44.000 align:middle line:84% Or you might have an entire chromosome, 00:07:44.000 --> 00:07:45.980 align:middle line:84% and you're asking questions about haplotype, 00:07:45.980 --> 00:07:48.170 align:middle line:84% as we were in the last class, which would be, 00:07:48.170 --> 00:07:54.854 align:middle line:84% what mutations are in cis on one chromosome relative to another? 00:07:54.854 --> 00:07:58.730 align:middle line:84% Another example is you might be scanning 00:07:58.730 --> 00:08:01.490 align:middle line:84% a chromosome for a motif that occurs again 00:08:01.490 --> 00:08:04.370 align:middle line:84% and again, like an Alu repeat or a transcription factor binding 00:08:04.370 --> 00:08:05.100 align:middle line:90% motif. 00:08:05.100 --> 00:08:07.100 align:middle line:84% And this might be in the middle of one sequence, 00:08:07.100 --> 00:08:09.740 align:middle line:84% at the end of another, or in the middle of both. 00:08:09.740 --> 00:08:12.180 align:middle line:84% Here, it's at the end of one and the middle of another. 00:08:12.180 --> 00:08:13.790 align:middle line:90% And it's short and local. 00:08:13.790 --> 00:08:16.070 align:middle line:84% And so you don't want to constrain 00:08:16.070 --> 00:08:18.050 align:middle line:90% the sequence ends to align. 00:08:18.050 --> 00:08:19.940 align:middle line:84% And you don't want to penalize if there's 00:08:19.940 --> 00:08:23.390 align:middle line:84% some deviation from alignment of the ends. 00:08:23.390 --> 00:08:27.710 align:middle line:84% Taking this one step further so that it's not internal 00:08:27.710 --> 00:08:32.630 align:middle line:84% to either one, in an ideal assembly-- 00:08:32.630 --> 00:08:35.809 align:middle line:84% shotgun assembly of sequence, as we talked about last time-- 00:08:35.809 --> 00:08:38.270 align:middle line:84% you would expect, as one fragment-- 00:08:38.270 --> 00:08:40.277 align:middle line:84% one clone fragment-- ends, hopefully, you 00:08:40.277 --> 00:08:42.860 align:middle line:84% have a little bit of overlap at the beginning of the next one, 00:08:42.860 --> 00:08:47.300 align:middle line:84% and you'll be able to jump along these stepping stones to get 00:08:47.300 --> 00:08:48.440 align:middle line:90% to the end of the sequence. 00:08:48.440 --> 00:08:52.250 align:middle line:84% Because you can't sequence the whole thing. 00:08:52.250 --> 00:08:55.700 align:middle line:84% But this ideal suffix, where a suffix of one 00:08:55.700 --> 00:09:03.365 align:middle line:84% overlaps the suffix, or the verse complement, of another, 00:09:03.365 --> 00:09:07.460 align:middle line:84% this sort of alignment can be imperfect 00:09:07.460 --> 00:09:10.080 align:middle line:84% because of errors at the ends of these sequences. 00:09:10.080 --> 00:09:13.010 align:middle line:84% So generally speaking, we'll be talking 00:09:13.010 --> 00:09:14.510 align:middle line:90% about global versus local. 00:09:14.510 --> 00:09:16.430 align:middle line:84% And these specific sequences we'll 00:09:16.430 --> 00:09:19.270 align:middle line:84% come back to several slides from now. 00:09:19.270 --> 00:09:21.020 align:middle line:84% Now, you want to have a scoring algorithm. 00:09:21.020 --> 00:09:27.140 align:middle line:84% We'll use this very simple scoring metric 00:09:27.140 --> 00:09:29.930 align:middle line:90% off and on during this talk. 00:09:29.930 --> 00:09:32.240 align:middle line:84% And here, we're just giving plus 1 00:09:32.240 --> 00:09:34.520 align:middle line:84% for every perfect match, indicated by a colon 00:09:34.520 --> 00:09:38.210 align:middle line:84% where an a matches an a, and a minus 1 for every mismatch, 00:09:38.210 --> 00:09:42.270 align:middle line:84% indicated by an x, where c does not match an a, for example. 00:09:42.270 --> 00:09:44.420 align:middle line:84% So we have five perfect matches and five 00:09:44.420 --> 00:09:46.080 align:middle line:90% mismatches in this case. 00:09:46.080 --> 00:09:51.440 align:middle line:84% And so it's a total of 5 minus 3 is 2. 00:09:51.440 --> 00:09:53.330 align:middle line:84% And in the local case, we're not going 00:09:53.330 --> 00:09:54.950 align:middle line:90% to require the ends to align. 00:09:54.950 --> 00:09:56.820 align:middle line:84% And so we can slip this a little bit 00:09:56.820 --> 00:09:59.180 align:middle line:84% so that now, there's four perfect matches, 00:09:59.180 --> 00:10:02.210 align:middle line:84% and we're not penalizing the terminal mismatches. 00:10:02.210 --> 00:10:05.030 align:middle line:84% And so we have a total of four, which, 00:10:05.030 --> 00:10:07.400 align:middle line:84% if these were directly comparable scores, 00:10:07.400 --> 00:10:10.010 align:middle line:90% would be a superior score. 00:10:10.010 --> 00:10:12.770 align:middle line:84% And you can see the suffixes, here, 00:10:12.770 --> 00:10:14.510 align:middle line:84% you're enforcing that it be at the ends. 00:10:14.510 --> 00:10:17.480 align:middle line:90% And it's a score of 3. 00:10:17.480 --> 00:10:20.500 align:middle line:84% So now we're going to compare different ways 00:10:20.500 --> 00:10:25.940 align:middle line:90% of searching through sequences. 00:10:25.940 --> 00:10:35.120 align:middle line:84% Now, that one was an exact match with mismatches, where 00:10:35.120 --> 00:10:36.590 align:middle line:90% the mismatches are penalized. 00:10:36.590 --> 00:10:41.480 align:middle line:84% And so an exact search, a truly exact search, is fairly rare-- 00:10:41.480 --> 00:10:43.820 align:middle line:84% maybe a restriction site, something like that. 00:10:43.820 --> 00:10:45.950 align:middle line:84% You can expand this to a regular expression 00:10:45.950 --> 00:10:48.770 align:middle line:84% where the insertions are restricted. 00:10:48.770 --> 00:10:53.450 align:middle line:84% They, in this case, an insertion can be indicated by any base. 00:10:53.450 --> 00:10:56.040 align:middle line:90% A, C, G, or T can occur. 00:10:56.040 --> 00:10:57.860 align:middle line:84% And then the number that you will tolerate 00:10:57.860 --> 00:11:00.370 align:middle line:84% is indicated by a numeric range-- zero 00:11:00.370 --> 00:11:01.620 align:middle line:90% to nine, in this case. 00:11:01.620 --> 00:11:04.670 align:middle line:84% And so the particular example given-- 00:11:04.670 --> 00:11:07.470 align:middle line:84% this equals sign just means that this is an example. 00:11:07.470 --> 00:11:09.920 align:middle line:90% So C, G are strict at the ends. 00:11:09.920 --> 00:11:13.340 align:middle line:84% And then it happened to have two As there, which 00:11:13.340 --> 00:11:15.770 align:middle line:84% an A is an example of a nucleotide that 00:11:15.770 --> 00:11:19.190 align:middle line:84% satisfies the abbreviation N. N is all possible. 00:11:19.190 --> 00:11:20.960 align:middle line:84% But you could similarly have the zero 00:11:20.960 --> 00:11:24.440 align:middle line:84% to nine bracket could refer to a short sequence, 00:11:24.440 --> 00:11:28.140 align:middle line:90% like a G sequence and so on. 00:11:28.140 --> 00:11:31.610 align:middle line:84% So you could get a known number of repeats. 00:11:31.610 --> 00:11:33.860 align:middle line:84% And that could represent the empirical observation 00:11:33.860 --> 00:11:35.990 align:middle line:84% that you have zero to nine AG repeats. 00:11:35.990 --> 00:11:40.920 align:middle line:90% 00:11:40.920 --> 00:11:43.500 align:middle line:84% Now, in the previous example, we just 00:11:43.500 --> 00:11:45.990 align:middle line:84% penalized whenever it was a mismatch by a fixed amount. 00:11:45.990 --> 00:11:47.920 align:middle line:84% But you could have a substitution matrix, 00:11:47.920 --> 00:11:51.870 align:middle line:84% which would actually codify your observations. 00:11:51.870 --> 00:11:54.700 align:middle line:84% When you look through natural sequences and you line them up, 00:11:54.700 --> 00:11:59.790 align:middle line:84% how often do you get an A substituting for a G? 00:11:59.790 --> 00:12:02.550 align:middle line:84% And so the penalty wouldn't be a strict minus 1. 00:12:02.550 --> 00:12:06.340 align:middle line:84% It would be something you'd look up in a table. 00:12:06.340 --> 00:12:10.560 align:middle line:84% And in that table, the diagonal would be the matches, 00:12:10.560 --> 00:12:13.200 align:middle line:84% and the off-diagonals would be specific mismatches. 00:12:13.200 --> 00:12:19.830 align:middle line:84% And we'll show how we get such a matrix and use it. 00:12:19.830 --> 00:12:22.740 align:middle line:84% Now, then you can have a profile matrix. 00:12:22.740 --> 00:12:27.720 align:middle line:84% Or this PSI means position-sensitive. 00:12:27.720 --> 00:12:31.140 align:middle line:84% And the position-sensitive means you have a different lookup 00:12:31.140 --> 00:12:34.320 align:middle line:84% table for every position in your biopolymer. 00:12:34.320 --> 00:12:37.080 align:middle line:84% And that makes sense because not all positions 00:12:37.080 --> 00:12:42.690 align:middle line:90% have the same set of equally-- 00:12:42.690 --> 00:12:46.090 align:middle line:84% or have the same sort of substitutions that are allowed. 00:12:46.090 --> 00:12:50.380 align:middle line:84% So in parentheses along here, these 00:12:50.380 --> 00:12:54.240 align:middle line:84% are actual names of programs that 00:12:54.240 --> 00:12:58.170 align:middle line:84% are available, either in commercial packages 00:12:58.170 --> 00:13:00.030 align:middle line:90% or for free. 00:13:00.030 --> 00:13:05.340 align:middle line:84% BLAST is Basic Local Alignment Sequence. 00:13:05.340 --> 00:13:08.850 align:middle line:90% And the N refers to nucleotide. 00:13:08.850 --> 00:13:11.795 align:middle line:84% And here, PSI-BLAST, again, position-sensitive. 00:13:11.795 --> 00:13:13.170 align:middle line:84% These are some of the ones you'll 00:13:13.170 --> 00:13:14.790 align:middle line:90% run into most frequently. 00:13:14.790 --> 00:13:19.650 align:middle line:90% 00:13:19.650 --> 00:13:23.640 align:middle line:84% The original versions of BLAST were basically 00:13:23.640 --> 00:13:30.120 align:middle line:84% aimed at having the largest block of contiguous sequence 00:13:30.120 --> 00:13:32.610 align:middle line:90% without gaps. 00:13:32.610 --> 00:13:35.910 align:middle line:84% And there were so-called Carlin statistics 00:13:35.910 --> 00:13:39.930 align:middle line:84% that would tell you that that's the largest sequence that will 00:13:39.930 --> 00:13:42.300 align:middle line:90% give you the best probability. 00:13:42.300 --> 00:13:44.910 align:middle line:84% But then it was widely recognized 00:13:44.910 --> 00:13:47.820 align:middle line:84% that when people actually evaluate matches between two 00:13:47.820 --> 00:13:51.390 align:middle line:84% sequences, they're not just evaluating the longest 00:13:51.390 --> 00:13:52.590 align:middle line:90% ungapped sequences. 00:13:52.590 --> 00:13:54.660 align:middle line:84% They're actually interested in the significance 00:13:54.660 --> 00:13:58.110 align:middle line:84% of a sequence that can include insertions and deletions. 00:13:58.110 --> 00:14:00.660 align:middle line:84% And so BLAST has been extended to include gaps. 00:14:00.660 --> 00:14:03.665 align:middle line:84% And prior to that long history back to the '70s 00:14:03.665 --> 00:14:05.415 align:middle line:84% is the Needleman-Wunsch and Smith-Waterman 00:14:05.415 --> 00:14:08.340 align:middle line:84% that I mentioned, a global and local, respectively, which 00:14:08.340 --> 00:14:13.560 align:middle line:84% would allow a large number of insertions and deletions 00:14:13.560 --> 00:14:15.660 align:middle line:84% of single bases and multiple bases, 00:14:15.660 --> 00:14:19.590 align:middle line:84% as determined by the parameters that were either manually 00:14:19.590 --> 00:14:21.960 align:middle line:90% set or determined empirically. 00:14:21.960 --> 00:14:23.610 align:middle line:84% And we'll try to stress ways that this 00:14:23.610 --> 00:14:25.350 align:middle line:90% can be determined empirically. 00:14:25.350 --> 00:14:29.850 align:middle line:84% And finally, the hidden Markov models truly 00:14:29.850 --> 00:14:35.760 align:middle line:84% take a probabilistic approach to each of these, 00:14:35.760 --> 00:14:43.440 align:middle line:84% allowing position-sensitive and models that are extendable 00:14:43.440 --> 00:14:48.510 align:middle line:84% not just where each position has a set of probabilities 00:14:48.510 --> 00:14:53.970 align:middle line:84% but there can be dependencies upon adjacent positions. 00:14:53.970 --> 00:14:57.870 align:middle line:84% And we'll get to this at the very end. 00:14:57.870 --> 00:15:00.180 align:middle line:84% And when we talked in the very first class 00:15:00.180 --> 00:15:02.820 align:middle line:84% about different definitions of complexity, 00:15:02.820 --> 00:15:05.670 align:middle line:84% one of them we talked about was the computational complexity, 00:15:05.670 --> 00:15:10.547 align:middle line:84% or hardness, of the amount of time that it takes, 00:15:10.547 --> 00:15:12.630 align:middle line:84% or the amount of space, or space and time, that it 00:15:12.630 --> 00:15:14.670 align:middle line:90% takes to solve a problem. 00:15:14.670 --> 00:15:19.200 align:middle line:84% And this certainly is a good illustration in today's talk 00:15:19.200 --> 00:15:21.750 align:middle line:90% in dynamic programming. 00:15:21.750 --> 00:15:25.763 align:middle line:84% Because when we want to do either a pairwise sequence 00:15:25.763 --> 00:15:27.430 align:middle line:84% alignment or a multisequence alignment-- 00:15:27.430 --> 00:15:29.280 align:middle line:84% let's start with pairwise, we're aligning k 00:15:29.280 --> 00:15:33.060 align:middle line:84% equals 2 sequences of length n, and we're allowing gaps. 00:15:33.060 --> 00:15:36.450 align:middle line:84% Now, if we're just comparing them without gaps, 00:15:36.450 --> 00:15:39.300 align:middle line:84% as we did in the earlier slides, it's trivial. 00:15:39.300 --> 00:15:42.288 align:middle line:84% It's linear with the total number of bases. 00:15:42.288 --> 00:15:43.080 align:middle line:90% It's linear within. 00:15:43.080 --> 00:15:45.000 align:middle line:84% And that, of course, scales very gracefully. 00:15:45.000 --> 00:15:48.360 align:middle line:84% But if you have a very naive-- and I'm going to set up a straw 00:15:48.360 --> 00:15:49.830 align:middle line:90% man here-- 00:15:49.830 --> 00:15:52.440 align:middle line:84% you have a very naive algorithm, then you'll go along, 00:15:52.440 --> 00:15:55.410 align:middle line:84% and you'll put in a gap at every possible position 00:15:55.410 --> 00:15:58.800 align:middle line:84% at the top in combination with every possible position 00:15:58.800 --> 00:15:59.770 align:middle line:90% in the bottom. 00:15:59.770 --> 00:16:01.670 align:middle line:84% And so that means that there's n on the top 00:16:01.670 --> 00:16:04.740 align:middle line:84% and n on the bottom for a total of 2n possible positions 00:16:04.740 --> 00:16:05.970 align:middle line:90% that you could put in a gap. 00:16:05.970 --> 00:16:08.170 align:middle line:84% And the gap can be any length, as you 00:16:08.170 --> 00:16:11.340 align:middle line:84% see on the right here of slide 7 is that you 00:16:11.340 --> 00:16:14.460 align:middle line:90% can have a gap up to length n. 00:16:14.460 --> 00:16:17.760 align:middle line:84% So roughly speaking, such a naive algorithm 00:16:17.760 --> 00:16:21.630 align:middle line:84% would scale very exponentially n to the 2n power. 00:16:21.630 --> 00:16:26.160 align:middle line:84% And this is just enormous for n the dozens. 00:16:26.160 --> 00:16:29.010 align:middle line:84% And when n is in the billions, then you can just basically 00:16:29.010 --> 00:16:31.140 align:middle line:90% forget about it. 00:16:31.140 --> 00:16:33.990 align:middle line:84% We're not even talking about computers or particles 00:16:33.990 --> 00:16:37.620 align:middle line:90% in the universe anymore. 00:16:37.620 --> 00:16:40.290 align:middle line:84% And that's for k equals 2, the simplest possible case where 00:16:40.290 --> 00:16:42.775 align:middle line:90% you're aligning two sequences. 00:16:42.775 --> 00:16:44.400 align:middle line:84% If you're aligning multiple sequences-- 00:16:44.400 --> 00:16:46.710 align:middle line:90% we'll get to that in a moment-- 00:16:46.710 --> 00:16:49.320 align:middle line:84% it definitely becomes exponential in k, even 00:16:49.320 --> 00:16:52.710 align:middle line:90% in a-non naive algorithm. 00:16:52.710 --> 00:16:54.600 align:middle line:84% So we're going to show that it's going 00:16:54.600 --> 00:16:58.720 align:middle line:84% to be non-exponential in n, the length of the sequences. 00:16:58.720 --> 00:17:04.380 align:middle line:84% But in order to do that, in order to assess algorithms, 00:17:04.380 --> 00:17:05.478 align:middle line:90% we can do them-- 00:17:05.478 --> 00:17:07.020 align:middle line:84% some of them we can do theoretically. 00:17:07.020 --> 00:17:10.000 align:middle line:84% And in fact, for the dynamic programming, 00:17:10.000 --> 00:17:12.030 align:middle line:84% we will show that you can convince yourself 00:17:12.030 --> 00:17:15.690 align:middle line:84% that an n squared algorithm is precise. 00:17:15.690 --> 00:17:18.030 align:middle line:84% But others, we will have to do empirically. 00:17:18.030 --> 00:17:19.619 align:middle line:84% And I just want to take an aside here 00:17:19.619 --> 00:17:27.079 align:middle line:84% to talk about how a particular comparison of sequence 00:17:27.079 --> 00:17:30.840 align:middle line:90% alignment programs was done. 00:17:30.840 --> 00:17:36.450 align:middle line:84% So the critical thing, not just for this test but for many that 00:17:36.450 --> 00:17:38.310 align:middle line:84% you'll be seeing in this course and you 00:17:38.310 --> 00:17:40.530 align:middle line:84% may want to do in your own projects, 00:17:40.530 --> 00:17:45.270 align:middle line:84% you want to set aside a training set in which you 00:17:45.270 --> 00:17:46.710 align:middle line:84% run through a number of algorithms 00:17:46.710 --> 00:17:48.660 align:middle line:84% or a number of parameters within an algorithm 00:17:48.660 --> 00:17:51.550 align:middle line:84% to assess which ones are the best. 00:17:51.550 --> 00:17:53.520 align:middle line:84% And once you've determined the algorithms that 00:17:53.520 --> 00:17:56.368 align:middle line:84% are best or the parameters that are best, 00:17:56.368 --> 00:17:58.660 align:middle line:84% then you want to have a testing set that's independent. 00:17:58.660 --> 00:18:02.070 align:middle line:84% That means non-redundant, as well. 00:18:02.070 --> 00:18:06.390 align:middle line:84% So the training set, if you use the training set 00:18:06.390 --> 00:18:08.970 align:middle line:84% as your testing set, you may lull yourself 00:18:08.970 --> 00:18:11.730 align:middle line:84% into a sense of complacency because it's been overtrained. 00:18:11.730 --> 00:18:14.520 align:middle line:84% And you're basically only capable of solving the problems 00:18:14.520 --> 00:18:15.690 align:middle line:90% that you set before it. 00:18:15.690 --> 00:18:17.400 align:middle line:84% So you want a separate testing set. 00:18:17.400 --> 00:18:19.630 align:middle line:90% And that was done in this case. 00:18:19.630 --> 00:18:22.780 align:middle line:84% We need some sort of evaluation criteria. 00:18:22.780 --> 00:18:25.410 align:middle line:84% So we talked a little bit of scores 00:18:25.410 --> 00:18:27.210 align:middle line:84% that you'd set up in order to score 00:18:27.210 --> 00:18:30.930 align:middle line:84% whether one alignment is better than another alignment. 00:18:30.930 --> 00:18:36.340 align:middle line:84% But in addition, when you want to compare two algorithms that 00:18:36.340 --> 00:18:41.130 align:middle line:84% may give you the same score, you want 00:18:41.130 --> 00:18:45.510 align:middle line:84% to have an external evaluation criterion. 00:18:45.510 --> 00:18:47.460 align:middle line:84% Typically, the evaluation criterion 00:18:47.460 --> 00:18:50.665 align:middle line:84% that one might have is sensitivity and specificity. 00:18:50.665 --> 00:18:52.290 align:middle line:84% Occasionally, you'll find in literature 00:18:52.290 --> 00:18:54.690 align:middle line:84% where people will just focus on one or the other. 00:18:54.690 --> 00:18:58.913 align:middle line:84% For some reason or other, they don't want to miss anything. 00:18:58.913 --> 00:19:00.330 align:middle line:84% So they want to keep their false-- 00:19:00.330 --> 00:19:02.970 align:middle line:84% they want to reduce their false positives as much as possible. 00:19:02.970 --> 00:19:06.640 align:middle line:84% Or they don't want to plow through mountains of output, 00:19:06.640 --> 00:19:08.940 align:middle line:84% and so they keep the false positives down. 00:19:08.940 --> 00:19:11.640 align:middle line:84% But you really want to have them both very low. 00:19:11.640 --> 00:19:14.760 align:middle line:84% And sometimes this is restated as sensitivity and specificity, 00:19:14.760 --> 00:19:23.140 align:middle line:84% where sensitivity is the number of true hits 00:19:23.140 --> 00:19:26.500 align:middle line:84% that you've predicted over the total number of true hits. 00:19:26.500 --> 00:19:28.450 align:middle line:84% This is, say, in your training set, 00:19:28.450 --> 00:19:30.220 align:middle line:84% where you know what the true hits are 00:19:30.220 --> 00:19:33.820 align:middle line:90% from some outside source. 00:19:33.820 --> 00:19:36.310 align:middle line:84% And then specificity is the number of true hits 00:19:36.310 --> 00:19:39.970 align:middle line:84% that you predicted over all of the hits that you predicted, 00:19:39.970 --> 00:19:41.870 align:middle line:90% true and false. 00:19:41.870 --> 00:19:45.080 align:middle line:84% Now, this truth here that comes in your training 00:19:45.080 --> 00:19:47.000 align:middle line:90% set, where does that come from? 00:19:47.000 --> 00:19:50.330 align:middle line:84% If we had access to truth, in general, then 00:19:50.330 --> 00:19:52.830 align:middle line:84% what do we need these algorithms for? 00:19:52.830 --> 00:19:55.760 align:middle line:84% And the answer is that we do have access 00:19:55.760 --> 00:19:59.690 align:middle line:84% to, maybe, a higher truth or something that's 00:19:59.690 --> 00:20:02.120 align:middle line:90% outside of sequence alignment. 00:20:02.120 --> 00:20:04.790 align:middle line:84% And this is crystallography, genetics, and biochemistry. 00:20:04.790 --> 00:20:06.290 align:middle line:84% Crystallography, for reasons we'll 00:20:06.290 --> 00:20:09.530 align:middle line:84% go into in the protein part of this course, 00:20:09.530 --> 00:20:12.410 align:middle line:84% is capable of detecting much more ancient relationships 00:20:12.410 --> 00:20:15.920 align:middle line:84% between biopolymers than is sequence alone. 00:20:15.920 --> 00:20:17.720 align:middle line:84% Similar genetics in biochemistry can 00:20:17.720 --> 00:20:23.120 align:middle line:84% test structural and functional hypotheses 00:20:23.120 --> 00:20:24.610 align:middle line:90% through great expense. 00:20:24.610 --> 00:20:25.790 align:middle line:90% And so these are expensive. 00:20:25.790 --> 00:20:28.580 align:middle line:84% And so that's the reason they're great for making training sets, 00:20:28.580 --> 00:20:31.010 align:middle line:84% but they won't necessarily replace 00:20:31.010 --> 00:20:36.380 align:middle line:90% sequence alignment and scores. 00:20:36.380 --> 00:20:38.630 align:middle line:84% So that was the setup for this slide, which 00:20:38.630 --> 00:20:42.070 align:middle line:84% is that Bill Pearson, who actually developed the FASTA 00:20:42.070 --> 00:20:49.880 align:middle line:84% algorithm, among others, did a thorough assessment 00:20:49.880 --> 00:20:51.260 align:middle line:90% of various algorithms. 00:20:51.260 --> 00:20:54.290 align:middle line:84% FASTA was one that was based on words, meaning 00:20:54.290 --> 00:20:58.610 align:middle line:84% exact matches of some fixed lengths the user could set; 00:20:58.610 --> 00:21:03.470 align:middle line:84% FASTP, which was based on these maximum-lengths blocks 00:21:03.470 --> 00:21:04.820 align:middle line:90% without gaps in it. 00:21:04.820 --> 00:21:08.587 align:middle line:84% Blitz is a variation on the Smith-Waterman algorithm, 00:21:08.587 --> 00:21:10.670 align:middle line:84% meaning a full dynamic programming, which we'll be 00:21:10.670 --> 00:21:12.980 align:middle line:90% talking about in a short while. 00:21:12.980 --> 00:21:14.750 align:middle line:84% This is a highly parallelized version 00:21:14.750 --> 00:21:17.010 align:middle line:84% of it, early parallel version of it. 00:21:17.010 --> 00:21:19.460 align:middle line:84% So the different algorithms for doing the alignment 00:21:19.460 --> 00:21:22.010 align:middle line:84% were compared, different substitution matrices 00:21:22.010 --> 00:21:23.030 align:middle line:90% and different databases. 00:21:23.030 --> 00:21:27.140 align:middle line:84% Just in case there was some database bias, 00:21:27.140 --> 00:21:30.240 align:middle line:90% he included that. 00:21:30.240 --> 00:21:32.830 align:middle line:84% And so we're going to talk about substitution matrices in just 00:21:32.830 --> 00:21:33.330 align:middle line:90% a moment. 00:21:33.330 --> 00:21:39.810 align:middle line:84% But basically, these are what amino acid or nucleotide-- 00:21:39.810 --> 00:21:42.180 align:middle line:84% what amino acid, in this case, can substitute 00:21:42.180 --> 00:21:47.880 align:middle line:84% for another amino acid in actual protein segments 00:21:47.880 --> 00:21:49.890 align:middle line:84% that have diverged by about the amount 00:21:49.890 --> 00:21:53.550 align:middle line:84% that you want to do your test on. 00:21:53.550 --> 00:21:55.440 align:middle line:84% And these different numbers just refer 00:21:55.440 --> 00:22:02.440 align:middle line:84% to roughly how distantly related the proteins are. 00:22:02.440 --> 00:22:04.500 align:middle line:84% The higher the number, the more distantly 00:22:04.500 --> 00:22:10.020 align:middle line:84% related they are in the case of the PAM matrices. 00:22:10.020 --> 00:22:14.480 align:middle line:84% So now, why did he do that with a protein level 00:22:14.480 --> 00:22:17.030 align:middle line:84% rather than doing it to the nucleic acid level? 00:22:17.030 --> 00:22:20.330 align:middle line:84% Well, historically, it's because there weren't any nucleic acid 00:22:20.330 --> 00:22:20.930 align:middle line:90% sequences. 00:22:20.930 --> 00:22:22.830 align:middle line:84% There were mostly protein sequences. 00:22:22.830 --> 00:22:29.330 align:middle line:84% But even today, there is obviously a lot more nucleic 00:22:29.330 --> 00:22:31.220 align:middle line:90% acid sequence. 00:22:31.220 --> 00:22:34.940 align:middle line:84% But there's a real reason to do it at the protein level, which 00:22:34.940 --> 00:22:37.100 align:middle line:84% is that when you look at the code 00:22:37.100 --> 00:22:41.720 align:middle line:84% that we've been talking about in these lectures, 00:22:41.720 --> 00:22:45.380 align:middle line:84% something like leucine can be represented 00:22:45.380 --> 00:22:48.560 align:middle line:84% by six different codons, which can have wildly 00:22:48.560 --> 00:22:50.270 align:middle line:90% different nucleotide sequences. 00:22:50.270 --> 00:22:56.180 align:middle line:84% So for example, CUG is valid, and UUA. 00:22:56.180 --> 00:23:00.980 align:middle line:84% And those only share one nucleotide out of three. 00:23:00.980 --> 00:23:03.320 align:middle line:84% And over long periods of time, if there's 00:23:03.320 --> 00:23:06.080 align:middle line:84% heavy selection on the protein and relatively weak selection 00:23:06.080 --> 00:23:08.480 align:middle line:84% on the nucleic acid, or there could even 00:23:08.480 --> 00:23:10.740 align:middle line:84% be pressure on the nucleic acid to change 00:23:10.740 --> 00:23:13.580 align:middle line:84% for reasons that we'll go into in the second half 00:23:13.580 --> 00:23:19.460 align:middle line:84% of this lecture, that pressure on the nucleotide sequence 00:23:19.460 --> 00:23:21.770 align:middle line:84% can cause the nucleotide sequence to change a lot 00:23:21.770 --> 00:23:23.930 align:middle line:84% and the protein sequence not to change much at all. 00:23:23.930 --> 00:23:26.880 align:middle line:90% 00:23:26.880 --> 00:23:30.113 align:middle line:84% So an example of pressure is if the tRNAs change 00:23:30.113 --> 00:23:32.030 align:middle line:84% in their abundance, then there'd be a pressure 00:23:32.030 --> 00:23:35.060 align:middle line:90% to codon usage to change. 00:23:35.060 --> 00:23:37.500 align:middle line:84% There are some reasons to do it at the nucleotide level. 00:23:37.500 --> 00:23:39.950 align:middle line:84% For example, if you're comparing sequences 00:23:39.950 --> 00:23:42.680 align:middle line:84% which don't encode proteins, that's an obvious reason. 00:23:42.680 --> 00:23:47.960 align:middle line:84% If you have a lot of insertion and deletions 00:23:47.960 --> 00:23:53.570 align:middle line:84% or a tricky biological phenomenon like RNA splicing 00:23:53.570 --> 00:23:57.170 align:middle line:84% that causes the protein to be out of phase-- 00:23:57.170 --> 00:24:00.740 align:middle line:90% 00:24:00.740 --> 00:24:04.040 align:middle line:84% the inferred protein to be out of phase or hard to infer-- 00:24:04.040 --> 00:24:06.500 align:middle line:84% then you would do it at the nucleotide sequence level. 00:24:06.500 --> 00:24:09.240 align:middle line:90% 00:24:09.240 --> 00:24:11.670 align:middle line:84% Now, I'm going to show this slide twice. 00:24:11.670 --> 00:24:14.460 align:middle line:84% The first time, we're going to take it as a given 00:24:14.460 --> 00:24:16.620 align:middle line:84% that we've been given this multisequence alignment, 00:24:16.620 --> 00:24:21.755 align:middle line:84% and we're not going to question right now how we got it. 00:24:21.755 --> 00:24:23.880 align:middle line:84% But we're going to use that multisequence alignment 00:24:23.880 --> 00:24:27.480 align:middle line:84% to derive, or to talk about how we would derive, 00:24:27.480 --> 00:24:29.980 align:middle line:90% a substitution matrix. 00:24:29.980 --> 00:24:32.373 align:middle line:84% And here, a substitution matrix, you 00:24:32.373 --> 00:24:34.290 align:middle line:84% can think-- this is a multisequence alignment. 00:24:34.290 --> 00:24:36.840 align:middle line:84% So essentially, we have a weight matrix, which, 00:24:36.840 --> 00:24:40.410 align:middle line:84% if this were position-sensitive, we would say, at this position, 00:24:40.410 --> 00:24:42.660 align:middle line:90% C never changes. 00:24:42.660 --> 00:24:44.400 align:middle line:84% If we do enough of these proteins 00:24:44.400 --> 00:24:46.470 align:middle line:84% and we don't care about position, 00:24:46.470 --> 00:24:48.850 align:middle line:84% we can build up a big set of matrices. 00:24:48.850 --> 00:24:51.210 align:middle line:84% And in general, we will find that C 00:24:51.210 --> 00:24:53.520 align:middle line:84% tends to substitute for C, and very few other things 00:24:53.520 --> 00:24:54.270 align:middle line:90% substitute for it. 00:24:54.270 --> 00:24:57.840 align:middle line:84% Eventually, you will find other substitutions. 00:24:57.840 --> 00:25:00.180 align:middle line:84% Cysteine and tryptophan, C and W, 00:25:00.180 --> 00:25:04.290 align:middle line:84% are relatively rare amino acids, and they're highly conserved. 00:25:04.290 --> 00:25:07.590 align:middle line:84% Other ones can be substituted, as you can see here-- 00:25:07.590 --> 00:25:12.010 align:middle line:84% threonine, serine, and valine can substitute for one another. 00:25:12.010 --> 00:25:14.520 align:middle line:84% So now let's take a look at how this plays out 00:25:14.520 --> 00:25:20.110 align:middle line:84% when we look at all the possible substitutions that can occur. 00:25:20.110 --> 00:25:22.470 align:middle line:90% And that's what's in slide 12. 00:25:22.470 --> 00:25:28.890 align:middle line:84% If you look along the very top row are the percentages, 00:25:28.890 --> 00:25:34.650 align:middle line:84% the abundances, of amino acids in a particular organism, 00:25:34.650 --> 00:25:36.510 align:middle line:90% say E. coli. 00:25:36.510 --> 00:25:43.920 align:middle line:84% And then there's a single-letter code, A through Y. 00:25:43.920 --> 00:25:48.870 align:middle line:84% And along the diagonal is the substitution matrix, 00:25:48.870 --> 00:25:50.670 align:middle line:84% the score which has been determined-- 00:25:50.670 --> 00:25:52.530 align:middle line:90% this is the BLOSUM matrix-- 00:25:52.530 --> 00:25:57.900 align:middle line:84% for a block represented in distantly related blocks. 00:25:57.900 --> 00:26:01.390 align:middle line:84% Here, you can see that the diagonal represents 00:26:01.390 --> 00:26:04.740 align:middle line:84% the tendency for the amino acid to substitute for itself. 00:26:04.740 --> 00:26:07.980 align:middle line:84% And those amino acids, which generally are not 00:26:07.980 --> 00:26:10.197 align:middle line:84% easily substituted-- as I say, are highly conserved, 00:26:10.197 --> 00:26:12.030 align:middle line:84% which we pointed out in the previous slide-- 00:26:12.030 --> 00:26:15.810 align:middle line:84% were cysteine, C, and tryptophan, W. 00:26:15.810 --> 00:26:19.230 align:middle line:84% And for example, W is the most strongly conserved. 00:26:19.230 --> 00:26:20.545 align:middle line:90% It's 22 along the diagonal. 00:26:20.545 --> 00:26:22.170 align:middle line:84% And the consequence of that is there'll 00:26:22.170 --> 00:26:26.860 align:middle line:84% be relatively few positive values off the diagonal. 00:26:26.860 --> 00:26:29.295 align:middle line:84% And in fact, for tryptophan, there are no positive values. 00:26:29.295 --> 00:26:31.980 align:middle line:90% 00:26:31.980 --> 00:26:35.610 align:middle line:84% The numbers here have been generated in such a way 00:26:35.610 --> 00:26:43.060 align:middle line:84% that the negatives will tend to cancel out 00:26:43.060 --> 00:26:48.580 align:middle line:84% the positives in alignments, known alignments, of sequences 00:26:48.580 --> 00:26:51.890 align:middle line:84% that are about the right evolutionary distance from one 00:26:51.890 --> 00:26:52.390 align:middle line:90% another. 00:26:52.390 --> 00:26:53.920 align:middle line:90% You want to pick-- 00:26:53.920 --> 00:26:56.470 align:middle line:84% you want to make your matrix, if you're 00:26:56.470 --> 00:26:58.540 align:middle line:84% trying to look for very distant related proteins, 00:26:58.540 --> 00:27:02.110 align:middle line:84% you want to take the substitutions 00:27:02.110 --> 00:27:04.360 align:middle line:84% that you're sampling in your training set 00:27:04.360 --> 00:27:09.250 align:middle line:84% to be as far at that same distance-- that 00:27:09.250 --> 00:27:10.850 align:middle line:90% is, very distantly related. 00:27:10.850 --> 00:27:13.550 align:middle line:84% And this is one of the mistakes that was made in the early PAM 00:27:13.550 --> 00:27:14.050 align:middle line:90% matrix. 00:27:14.050 --> 00:27:17.110 align:middle line:84% There were two mistakes, actually-- the mathematics-- 00:27:17.110 --> 00:27:20.050 align:middle line:84% well, first, the proteins that were compared 00:27:20.050 --> 00:27:21.190 align:middle line:90% were very closely related. 00:27:21.190 --> 00:27:23.650 align:middle line:84% Because closely related proteins were more trustworthy. 00:27:23.650 --> 00:27:26.920 align:middle line:84% You could align closed sequences more easily. 00:27:26.920 --> 00:27:29.800 align:middle line:84% The algorithms didn't have to be sophisticated. 00:27:29.800 --> 00:27:32.000 align:middle line:84% And the trees could be more precise. 00:27:32.000 --> 00:27:34.060 align:middle line:84% But then, because they had done-- 00:27:34.060 --> 00:27:37.720 align:middle line:84% so that already was a bias because the substitutions 00:27:37.720 --> 00:27:39.220 align:middle line:84% you get in closely related sequences 00:27:39.220 --> 00:27:40.300 align:middle line:90% aren't really the same. 00:27:40.300 --> 00:27:42.850 align:middle line:84% And then they apply to a mathematical extrapolation 00:27:42.850 --> 00:27:46.270 align:middle line:84% method, which was not adequate in terms 00:27:46.270 --> 00:27:49.900 align:middle line:84% of the actual evolution and also wasn't even correct, 00:27:49.900 --> 00:27:54.100 align:middle line:84% mathematically, although this persisted for decades 00:27:54.100 --> 00:27:54.970 align:middle line:90% as the most common-- 00:27:54.970 --> 00:27:57.010 align:middle line:84% and still-- the most commonly used matrix. 00:27:57.010 --> 00:27:59.980 align:middle line:90% 00:27:59.980 --> 00:28:02.410 align:middle line:84% Anyway, you can see that, although tryptophan doesn't 00:28:02.410 --> 00:28:04.840 align:middle line:84% have any positive off-diagonal, something like arginine, 00:28:04.840 --> 00:28:10.750 align:middle line:84% here, in this blue, has a positive 4 off-diagonal 00:28:10.750 --> 00:28:14.380 align:middle line:90% and a positive 10 on-diagonal. 00:28:14.380 --> 00:28:20.230 align:middle line:84% So as you might guess, what's the most likely 00:28:20.230 --> 00:28:23.050 align:middle line:84% substituted for positively charged arginine? 00:28:23.050 --> 00:28:27.172 align:middle line:84% It's positively charged lysine under physiological conditions. 00:28:27.172 --> 00:28:28.630 align:middle line:84% And that's why that's off-diagonal. 00:28:28.630 --> 00:28:29.672 align:middle line:90% And there are other ones. 00:28:29.672 --> 00:28:32.620 align:middle line:84% We've color coded these the same as the genetic code, where 00:28:32.620 --> 00:28:37.390 align:middle line:84% the negatively charged amino acids can also 00:28:37.390 --> 00:28:38.590 align:middle line:90% substitute for one another. 00:28:38.590 --> 00:28:41.740 align:middle line:90% 00:28:41.740 --> 00:28:46.120 align:middle line:84% So the significance of that top row, of the percent abundance, 00:28:46.120 --> 00:28:48.640 align:middle line:84% is that if you find two matching As, 00:28:48.640 --> 00:28:50.470 align:middle line:84% that's not so significant because that's 00:28:50.470 --> 00:28:54.430 align:middle line:84% the most frequently occurring amino acid in this organism. 00:28:54.430 --> 00:28:58.030 align:middle line:84% On the other hand, if you find two matching Cs, 00:28:58.030 --> 00:28:59.560 align:middle line:90% that's very significant. 00:28:59.560 --> 00:29:00.880 align:middle line:90% Because those are rare. 00:29:00.880 --> 00:29:03.250 align:middle line:84% And finding two of them at the same place 00:29:03.250 --> 00:29:06.800 align:middle line:84% in a particular alignment means it's significant. 00:29:06.800 --> 00:29:10.990 align:middle line:84% So both the abundance and the substitution matrix 00:29:10.990 --> 00:29:13.610 align:middle line:90% can be useful. 00:29:13.610 --> 00:29:19.090 align:middle line:84% So now, we're going to walk through a actual scoring 00:29:19.090 --> 00:29:21.190 align:middle line:90% of some alignments. 00:29:21.190 --> 00:29:26.260 align:middle line:84% And we want to do this in this more challenging situation 00:29:26.260 --> 00:29:28.750 align:middle line:84% where you allow insertions and deletions. 00:29:28.750 --> 00:29:33.610 align:middle line:84% So first, we will, even though we've 00:29:33.610 --> 00:29:41.470 align:middle line:84% told how it is that we get the match versus mismatch numbers 00:29:41.470 --> 00:29:44.560 align:middle line:84% as a full substitution matrix, here, you 00:29:44.560 --> 00:29:48.220 align:middle line:84% can imagine substitution matrix has plus 1's along its diagonal 00:29:48.220 --> 00:29:50.068 align:middle line:84% and minus 1's on off-diagonal just so you 00:29:50.068 --> 00:29:51.610 align:middle line:84% can do all the calculations I'm going 00:29:51.610 --> 00:29:54.100 align:middle line:84% to show you in the next few slides in your head. 00:29:54.100 --> 00:29:56.127 align:middle line:84% But also, simultaneously imagine that it 00:29:56.127 --> 00:29:58.210 align:middle line:84% could have the richness of the substitution matrix 00:29:58.210 --> 00:29:59.057 align:middle line:90% we just had. 00:29:59.057 --> 00:30:01.390 align:middle line:84% We're going to do this, the next few, with nucleic acid. 00:30:01.390 --> 00:30:04.390 align:middle line:84% But imagine you could also do it with amino acids. 00:30:04.390 --> 00:30:07.960 align:middle line:84% The nucleic acid substitution matrix will be a 4 by 4. 00:30:07.960 --> 00:30:12.220 align:middle line:84% The one we just saw was a 20 by 20 for amino acids. 00:30:12.220 --> 00:30:16.450 align:middle line:84% The indels, we'll penalize them by minus 2. 00:30:16.450 --> 00:30:18.820 align:middle line:90% But this is an arbitrary number. 00:30:18.820 --> 00:30:22.700 align:middle line:84% And you'll see how critical it is in just a moment. 00:30:22.700 --> 00:30:25.420 align:middle line:84% But you can imagine that this could be determined 00:30:25.420 --> 00:30:27.580 align:middle line:84% empirically, just like the substitution 00:30:27.580 --> 00:30:29.200 align:middle line:84% matrix was determined empirically 00:30:29.200 --> 00:30:32.290 align:middle line:90% in the previous slide. 00:30:32.290 --> 00:30:35.230 align:middle line:84% The alignment score, then, will be defined as a sum of columns. 00:30:35.230 --> 00:30:38.770 align:middle line:84% We're going to be assuming that adjacent positions are 00:30:38.770 --> 00:30:40.610 align:middle line:90% independent of one another. 00:30:40.610 --> 00:30:43.180 align:middle line:84% And we'll be scoring them independently and then 00:30:43.180 --> 00:30:44.650 align:middle line:90% just taking the sum. 00:30:44.650 --> 00:30:46.360 align:middle line:84% That gives us the alignment score 00:30:46.360 --> 00:30:51.760 align:middle line:84% for a particular alignment, a particular set of indels 00:30:51.760 --> 00:30:55.180 align:middle line:84% and a particular set of offsets from one 00:30:55.180 --> 00:30:56.902 align:middle line:90% sequence relative to another. 00:30:56.902 --> 00:30:58.360 align:middle line:84% But what we really want to do is go 00:30:58.360 --> 00:31:00.310 align:middle line:84% through all those possible alignments 00:31:00.310 --> 00:31:01.990 align:middle line:84% to get the optimal alignment, which is 00:31:01.990 --> 00:31:05.403 align:middle line:90% the maximal score defined here. 00:31:05.403 --> 00:31:06.820 align:middle line:84% To get the optimal alignment, we'd 00:31:06.820 --> 00:31:10.270 align:middle line:84% like to do that in less than exponential time with at least 00:31:10.270 --> 00:31:13.840 align:middle line:84% less than exponential per n length of the sequences. 00:31:13.840 --> 00:31:16.870 align:middle line:90% 00:31:16.870 --> 00:31:25.790 align:middle line:84% So we're going to use this pair of sequences, ATGA, ACTA, 00:31:25.790 --> 00:31:26.810 align:middle line:90% twice on this slide. 00:31:26.810 --> 00:31:28.727 align:middle line:84% And we're going to use it in subsequent slides 00:31:28.727 --> 00:31:31.100 align:middle line:84% where we do it a slightly different way. 00:31:31.100 --> 00:31:34.460 align:middle line:84% And we're going to use that very simple scoring metric-- 00:31:34.460 --> 00:31:38.480 align:middle line:84% plus 1 from that perfect match, minus 1 for mismatch, minus 2 00:31:38.480 --> 00:31:40.520 align:middle line:90% for indel. 00:31:40.520 --> 00:31:45.080 align:middle line:84% And what we get here is these are just 00:31:45.080 --> 00:31:46.760 align:middle line:84% two of many different alignments we 00:31:46.760 --> 00:31:49.310 align:middle line:84% could have with different insertions here. 00:31:49.310 --> 00:31:52.250 align:middle line:84% On the left, number one, the most extreme case, 00:31:52.250 --> 00:31:54.498 align:middle line:84% no insertions or deletions on either sequence. 00:31:54.498 --> 00:31:55.790 align:middle line:90% We're only counting mismatches. 00:31:55.790 --> 00:31:57.620 align:middle line:84% There are two matches, two mismatches. 00:31:57.620 --> 00:32:00.590 align:middle line:84% And so that's two plus 1's, two minus 1's. 00:32:00.590 --> 00:32:01.490 align:middle line:90% And they cancel out. 00:32:01.490 --> 00:32:05.640 align:middle line:84% And the score is 0 for the one on the left. 00:32:05.640 --> 00:32:09.080 align:middle line:84% Then the one on the right, we've gone to the-- 00:32:09.080 --> 00:32:17.190 align:middle line:84% allowed an insertion on each strand, indicated by a dash 00:32:17.190 --> 00:32:20.520 align:middle line:90% on the opposite sequence. 00:32:20.520 --> 00:32:24.640 align:middle line:84% And now, you see you have three perfect matches, which 00:32:24.640 --> 00:32:26.640 align:middle line:84% is an increase in the number of perfect matches, 00:32:26.640 --> 00:32:31.710 align:middle line:84% but penalized with two indels, which are both negative 2. 00:32:31.710 --> 00:32:35.250 align:middle line:84% So it's plus 2 minus 2 minus 2 for minus 1. 00:32:35.250 --> 00:32:41.370 align:middle line:84% So this is not an improvement for the alignment on the left 00:32:41.370 --> 00:32:45.600 align:middle line:84% if we accept the scoring metrics that we 00:32:45.600 --> 00:32:47.310 align:middle line:90% had in the previous slide-- 00:32:47.310 --> 00:32:49.590 align:middle line:84% plus 1, minus 1, and minus 2 for the indels. 00:32:49.590 --> 00:32:52.260 align:middle line:84% If, instead, we say, well, indels really shouldn't be 00:32:52.260 --> 00:32:55.440 align:middle line:84% penalized that much-- we can accept insertions and deletions 00:32:55.440 --> 00:32:58.860 align:middle line:84% in these kind of sequences, we'll penalize by minus 1, 00:32:58.860 --> 00:33:01.110 align:middle line:84% penalize it the same as a mismatch-- 00:33:01.110 --> 00:33:07.960 align:middle line:90% now the score is the 3-- 00:33:07.960 --> 00:33:11.200 align:middle line:84% score for the one on the right, with three perfect matches 00:33:11.200 --> 00:33:14.110 align:middle line:84% and two insertion-deletions, is now plus 1. 00:33:14.110 --> 00:33:17.470 align:middle line:84% And it beats the perfectly aligned ones. 00:33:17.470 --> 00:33:22.335 align:middle line:84% So whether 1 is better than 2 depends on the indel score 00:33:22.335 --> 00:33:22.960 align:middle line:90% that you chose. 00:33:22.960 --> 00:33:24.377 align:middle line:84% STUDENT: Can I ask you a question? 00:33:24.377 --> 00:33:25.780 align:middle line:90% GEORGE CHURCH: Yes, question. 00:33:25.780 --> 00:33:27.197 align:middle line:84% STUDENT: Now that you are aligning 00:33:27.197 --> 00:33:32.040 align:middle line:84% two different sequences and in the case of the indel, 00:33:32.040 --> 00:33:35.440 align:middle line:84% you are allowing insertions all over the place-- 00:33:35.440 --> 00:33:38.520 align:middle line:84% I mean, you could theoretically have millions of those. 00:33:38.520 --> 00:33:42.270 align:middle line:84% But in reality, this [INAUDIBLE] the sequences in most cases 00:33:42.270 --> 00:33:43.160 align:middle line:90% would be no. 00:33:43.160 --> 00:33:44.300 align:middle line:90% You would know what it is. 00:33:44.300 --> 00:33:44.630 align:middle line:90% You can't-- 00:33:44.630 --> 00:33:46.588 align:middle line:84% GEORGE CHURCH: We know both of these sequences. 00:33:46.588 --> 00:33:50.300 align:middle line:84% These are sequences that we're comparing 00:33:50.300 --> 00:33:51.485 align:middle line:90% two different organisms. 00:33:51.485 --> 00:33:52.260 align:middle line:90% STUDENT: Right. 00:33:52.260 --> 00:33:54.710 align:middle line:84% So if you know both of them, then what's 00:33:54.710 --> 00:33:58.520 align:middle line:84% the point in allowing all these indel [INAUDIBLE].. 00:33:58.520 --> 00:34:01.370 align:middle line:84% GEORGE CHURCH: So the question is, why are we allowing 00:34:01.370 --> 00:34:02.570 align:middle line:90% insertions and deletions? 00:34:02.570 --> 00:34:07.430 align:middle line:84% And the reason is that during evolution, say, either lab 00:34:07.430 --> 00:34:12.290 align:middle line:84% evolution or ancient, insertions and deletions 00:34:12.290 --> 00:34:13.980 align:middle line:90% are valid mutations. 00:34:13.980 --> 00:34:17.690 align:middle line:84% And so we're trying to determine where the most likely 00:34:17.690 --> 00:34:20.600 align:middle line:84% places that insertions and deletions might 00:34:20.600 --> 00:34:24.199 align:middle line:84% have occurred over the course of the divergence of two 00:34:24.199 --> 00:34:25.310 align:middle line:90% sequences. 00:34:25.310 --> 00:34:27.500 align:middle line:84% And believe me, insertions and deletions 00:34:27.500 --> 00:34:29.909 align:middle line:90% are very, very common. 00:34:29.909 --> 00:34:32.210 align:middle line:90% So that's why we permit them. 00:34:32.210 --> 00:34:36.290 align:middle line:84% Now, why it is that insertion-deletions 00:34:36.290 --> 00:34:39.290 align:middle line:84% might be highly penalized or low penalized 00:34:39.290 --> 00:34:41.699 align:middle line:84% might depend on a position in the sequence. 00:34:41.699 --> 00:34:45.860 align:middle line:84% So for example, if you have a transcription factor where 00:34:45.860 --> 00:34:49.340 align:middle line:84% its precise geometry is important 00:34:49.340 --> 00:34:54.650 align:middle line:84% or an alpha helix and a protein or the translation 00:34:54.650 --> 00:34:56.989 align:middle line:84% of a genetic code where an insertion will 00:34:56.989 --> 00:35:01.730 align:middle line:84% throw the entire frame out of whack, as we 00:35:01.730 --> 00:35:06.200 align:middle line:84% had in the chemokine receptor in last class, 00:35:06.200 --> 00:35:09.140 align:middle line:84% then you want to penalize an indel very heavily. 00:35:09.140 --> 00:35:10.850 align:middle line:84% On the other hand, if you have a bunch 00:35:10.850 --> 00:35:15.303 align:middle line:84% of motifs that are kind of separated by variable linkers, 00:35:15.303 --> 00:35:16.970 align:middle line:84% then the insertion-deletion could almost 00:35:16.970 --> 00:35:18.720 align:middle line:90% be zero, no penalty at all. 00:35:18.720 --> 00:35:21.620 align:middle line:84% So you can see it matters, and it might be position-sensitive. 00:35:21.620 --> 00:35:24.560 align:middle line:84% It might not be one size fits all. 00:35:24.560 --> 00:35:27.230 align:middle line:84% But these are empirically determined-- 00:35:27.230 --> 00:35:28.820 align:middle line:90% can be empirically determined. 00:35:28.820 --> 00:35:32.060 align:middle line:84% So here's the hero-- dynamic programming. 00:35:32.060 --> 00:35:35.450 align:middle line:84% We've hopefully motivated that we can do scoring. 00:35:35.450 --> 00:35:39.560 align:middle line:84% We can determine empirically useful substitution matrices 00:35:39.560 --> 00:35:40.850 align:middle line:90% and indels. 00:35:40.850 --> 00:35:43.670 align:middle line:90% Now, how do we apply them? 00:35:43.670 --> 00:35:49.700 align:middle line:84% And dynamic programming extends beyond biology, 00:35:49.700 --> 00:35:51.590 align:middle line:90% as I've alluded. 00:35:51.590 --> 00:35:56.450 align:middle line:84% And such an algorithm solves every subproblem just once 00:35:56.450 --> 00:35:59.420 align:middle line:84% and saves its answer in a table, thereby avoiding 00:35:59.420 --> 00:36:02.610 align:middle line:84% the work of recomputing the answer every time. 00:36:02.610 --> 00:36:06.440 align:middle line:84% So the strawman that I threw up before 00:36:06.440 --> 00:36:11.180 align:middle line:84% of having this exponential problem is very readily solved. 00:36:11.180 --> 00:36:14.720 align:middle line:84% And the way it's solved is this is the subproblem 00:36:14.720 --> 00:36:16.050 align:middle line:90% way of dealing with it. 00:36:16.050 --> 00:36:18.980 align:middle line:84% And the idea of recursion, which we lightly 00:36:18.980 --> 00:36:21.830 align:middle line:84% touched upon when we defined the factorial, 00:36:21.830 --> 00:36:26.130 align:middle line:84% n factorial as equal to n times n minus 1 factorial, 00:36:26.130 --> 00:36:28.640 align:middle line:84% so defining it in terms of itself. 00:36:28.640 --> 00:36:31.520 align:middle line:84% But the key thing behind that definition and the ones 00:36:31.520 --> 00:36:34.130 align:middle line:84% we'll have here is that when you define something 00:36:34.130 --> 00:36:36.230 align:middle line:84% it terms of itself, you'd better have 00:36:36.230 --> 00:36:41.600 align:middle line:84% the call be a simpler problem and eventually terminate. 00:36:41.600 --> 00:36:44.600 align:middle line:90% And so that's what we have here. 00:36:44.600 --> 00:36:48.320 align:middle line:84% I'm going to give two examples in slide 17 00:36:48.320 --> 00:36:52.670 align:middle line:84% and in the next slide, one of them global and one of them 00:36:52.670 --> 00:36:53.510 align:middle line:90% local. 00:36:53.510 --> 00:36:55.160 align:middle line:84% This one will be done in terms of 00:36:55.160 --> 00:36:57.470 align:middle line:84% a tetranucleotide comparison, the same one we've 00:36:57.470 --> 00:36:58.655 align:middle line:90% been dealing with all along. 00:36:58.655 --> 00:37:00.905 align:middle line:84% And the other one will be in a more abstract sequence. 00:37:00.905 --> 00:37:04.940 align:middle line:90% 00:37:04.940 --> 00:37:10.220 align:middle line:84% Here, the way do the sub-subproblem by recursion is 00:37:10.220 --> 00:37:14.930 align:middle line:84% we say we define the score of aligning these two 00:37:14.930 --> 00:37:19.220 align:middle line:84% tetranucleotides as the maximum of-- 00:37:19.220 --> 00:37:20.810 align:middle line:84% and then there are three options. 00:37:20.810 --> 00:37:24.650 align:middle line:84% It can be either the score of having 00:37:24.650 --> 00:37:32.060 align:middle line:84% an insertion on the top strand, and that's the top option. 00:37:32.060 --> 00:37:34.760 align:middle line:84% The middle one is having no insertions or deletions 00:37:34.760 --> 00:37:37.640 align:middle line:84% on either strand and just evaluating 00:37:37.640 --> 00:37:43.810 align:middle line:84% the last base comparison, which, in this case, 00:37:43.810 --> 00:37:45.280 align:middle line:90% is an A versus an A. 00:37:45.280 --> 00:37:48.400 align:middle line:84% Now, that is the way that the algorithm terminates. 00:37:48.400 --> 00:37:52.240 align:middle line:84% When you have a single-base comparison or a single base 00:37:52.240 --> 00:37:55.600 align:middle line:84% compared to an indel, then you look 00:37:55.600 --> 00:37:59.200 align:middle line:84% up the scoring algorithm, or the scoring metrics 00:37:59.200 --> 00:38:00.500 align:middle line:90% we've been using all along. 00:38:00.500 --> 00:38:04.090 align:middle line:84% So here, let's look at that final right-hand column. 00:38:04.090 --> 00:38:08.560 align:middle line:84% The score for an indel versus an A would be that minus 2 00:38:08.560 --> 00:38:10.060 align:middle line:90% that we've been assuming along. 00:38:10.060 --> 00:38:11.890 align:middle line:84% And the score of an A versus A would 00:38:11.890 --> 00:38:13.900 align:middle line:84% be the substitution matrix diagonal, 00:38:13.900 --> 00:38:17.140 align:middle line:84% which would be a plus 1 and then, here, a minus 2. 00:38:17.140 --> 00:38:20.770 align:middle line:84% And so you can see that you're calling up 00:38:20.770 --> 00:38:26.440 align:middle line:84% these three possibilities-- indel, no insertion on the top, 00:38:26.440 --> 00:38:28.595 align:middle line:90% no insertion on the bottom. 00:38:28.595 --> 00:38:30.970 align:middle line:84% And you take the maximum of these, whichever one of these 00:38:30.970 --> 00:38:32.980 align:middle line:90% gives the best score. 00:38:32.980 --> 00:38:35.470 align:middle line:84% Now, that requires going back and calling it again. 00:38:35.470 --> 00:38:38.590 align:middle line:84% But you're calling it with a simpler-- 00:38:38.590 --> 00:38:40.370 align:middle line:90% you're asking for a simpler one. 00:38:40.370 --> 00:38:43.820 align:middle line:84% So now, you'll take the max of ATG versus ACT. 00:38:43.820 --> 00:38:47.230 align:middle line:84% And that'll ask you to look up the max of AT AC. 00:38:47.230 --> 00:38:49.480 align:middle line:84% And finally, it'll get the max of A versus A. 00:38:49.480 --> 00:38:51.220 align:middle line:90% And then you end. 00:38:51.220 --> 00:38:52.277 align:middle line:90% STUDENT: Excuse me. 00:38:52.277 --> 00:38:53.110 align:middle line:90% GEORGE CHURCH: Yeah. 00:38:53.110 --> 00:38:56.310 align:middle line:84% STUDENT: Are you assuming that insertion side is always a 1? 00:38:56.310 --> 00:38:59.830 align:middle line:84% The insertion side is always 1, right? 00:38:59.830 --> 00:39:00.580 align:middle line:90% GEORGE CHURCH: No. 00:39:00.580 --> 00:39:02.890 align:middle line:84% This algorithm allows any number of insertions 00:39:02.890 --> 00:39:04.450 align:middle line:90% up to the length of a sequence. 00:39:04.450 --> 00:39:06.790 align:middle line:84% And you'll see it when we do this in tabular form, 00:39:06.790 --> 00:39:09.250 align:middle line:90% how every possibility. 00:39:09.250 --> 00:39:11.000 align:middle line:90% But you do one at a time. 00:39:11.000 --> 00:39:12.790 align:middle line:90% There are only three cases here. 00:39:12.790 --> 00:39:15.250 align:middle line:84% By dealing with just three cases at a time, 00:39:15.250 --> 00:39:19.163 align:middle line:84% you actually end up having the full generality of any number 00:39:19.163 --> 00:39:20.330 align:middle line:90% of insertions and deletions. 00:39:20.330 --> 00:39:21.997 align:middle line:84% And that's the beauty of this algorithm. 00:39:21.997 --> 00:39:25.360 align:middle line:84% You don't have to explicitly do every possible insertion 00:39:25.360 --> 00:39:27.520 align:middle line:90% with every possible deletion. 00:39:27.520 --> 00:39:30.340 align:middle line:84% You just have to run them through once. 00:39:30.340 --> 00:39:31.890 align:middle line:90% OK. 00:39:31.890 --> 00:39:35.140 align:middle line:84% Now, I said that I was going to do two treatments 00:39:35.140 --> 00:39:38.170 align:middle line:90% of certain similarities. 00:39:38.170 --> 00:39:40.210 align:middle line:84% These are both dynamic programming 00:39:40.210 --> 00:39:42.910 align:middle line:90% of pairwise sequences. 00:39:42.910 --> 00:39:44.560 align:middle line:90% The previous one was global. 00:39:44.560 --> 00:39:45.610 align:middle line:90% This one is local. 00:39:45.610 --> 00:39:52.090 align:middle line:84% The only difference now is that we restrain the score 00:39:52.090 --> 00:39:53.470 align:middle line:90% to be greater than 0. 00:39:53.470 --> 00:39:55.190 align:middle line:90% We don't permit negative. 00:39:55.190 --> 00:39:57.370 align:middle line:84% So that means we're not penalizing the mismatches, 00:39:57.370 --> 00:39:58.510 align:middle line:90% for example, at the ends. 00:39:58.510 --> 00:40:01.340 align:middle line:84% Remember when I showed you that specific example early on. 00:40:01.340 --> 00:40:03.040 align:middle line:90% So now we have four choices-- 00:40:03.040 --> 00:40:07.258 align:middle line:90% the same three as before plus 0. 00:40:07.258 --> 00:40:09.550 align:middle line:84% And the other thing that I made a little different here 00:40:09.550 --> 00:40:12.550 align:middle line:84% is, rather than having a specific sequence-- that tetra 00:40:12.550 --> 00:40:14.830 align:middle line:84% nucleotide-- here, we have a general sequence where 00:40:14.830 --> 00:40:20.290 align:middle line:84% you show that the ellipses, the sequence is up to i 00:40:20.290 --> 00:40:22.810 align:middle line:90% long and up to j long here. 00:40:22.810 --> 00:40:25.990 align:middle line:84% And at this stage in the scoring, 00:40:25.990 --> 00:40:34.270 align:middle line:84% you're going to either lop off the i and j sequence element-- 00:40:34.270 --> 00:40:37.330 align:middle line:84% this would be a single base, and you do that score 00:40:37.330 --> 00:40:39.490 align:middle line:90% in the central scoring here-- 00:40:39.490 --> 00:40:45.530 align:middle line:84% or you have an insertion at the top or insertion at the bottom. 00:40:45.530 --> 00:40:50.410 align:middle line:84% So this is just a restatement of the previous one, generalized 00:40:50.410 --> 00:40:52.990 align:middle line:84% and made into a local alignment, which, 00:40:52.990 --> 00:40:54.820 align:middle line:90% in general, is what people do. 00:40:54.820 --> 00:40:57.790 align:middle line:84% People do local alignments rather than global ones. 00:40:57.790 --> 00:41:00.910 align:middle line:84% Because it's unsafe to say that the ends of your sequence 00:41:00.910 --> 00:41:01.600 align:middle line:90% will align. 00:41:01.600 --> 00:41:06.430 align:middle line:84% But we'll work through both of these as examples. 00:41:06.430 --> 00:41:13.860 align:middle line:84% Now, we're going to compute this as a row-by-row algorithm. 00:41:13.860 --> 00:41:18.630 align:middle line:84% Now, casually, you could just leave off 00:41:18.630 --> 00:41:22.930 align:middle line:90% this frame along the edge. 00:41:22.930 --> 00:41:27.390 align:middle line:84% But in order to make the algorithm be 00:41:27.390 --> 00:41:30.480 align:middle line:84% the same for the beginning and all the intermediate steps, 00:41:30.480 --> 00:41:34.950 align:middle line:84% what you do is you pre-fill this with numbers 00:41:34.950 --> 00:41:43.320 align:middle line:84% such that the edges are some very, very small number which 00:41:43.320 --> 00:41:45.670 align:middle line:84% is smaller than the sum of all the scores 00:41:45.670 --> 00:41:47.170 align:middle line:84% that you could get out of this table 00:41:47.170 --> 00:41:50.970 align:middle line:84% so that you can't really come in from those edges. 00:41:50.970 --> 00:41:54.030 align:middle line:84% You have to come in from the zero 00:41:54.030 --> 00:41:57.060 align:middle line:84% point, which is the global alignment requires 00:41:57.060 --> 00:41:57.960 align:middle line:90% the ends to align. 00:41:57.960 --> 00:42:00.900 align:middle line:84% So this is requiring the left-hand end to align. 00:42:00.900 --> 00:42:03.445 align:middle line:84% And so then the first comparison-- 00:42:03.445 --> 00:42:05.070 align:middle line:84% the only comparison you can really do-- 00:42:05.070 --> 00:42:07.715 align:middle line:84% is the A, A. That's the terminal comparison. 00:42:07.715 --> 00:42:09.340 align:middle line:84% And that happens to be a perfect match, 00:42:09.340 --> 00:42:10.610 align:middle line:90% so it gets a score of plus 1. 00:42:10.610 --> 00:42:15.960 align:middle line:84% Now, the next square that you can do is minus 1. 00:42:15.960 --> 00:42:21.330 align:middle line:84% And remember, each of these has three possibilities 00:42:21.330 --> 00:42:22.330 align:middle line:90% in the global alignment. 00:42:22.330 --> 00:42:27.260 align:middle line:84% Remember, it can be an insertion, a deletion, 00:42:27.260 --> 00:42:31.550 align:middle line:84% or just a direct comparison of match versus mismatch. 00:42:31.550 --> 00:42:34.760 align:middle line:84% So this first one, the insertion and deletion were ruled out. 00:42:34.760 --> 00:42:38.990 align:middle line:84% They weren't going to win the maximum score. 00:42:38.990 --> 00:42:41.065 align:middle line:90% So you basically got 1. 00:42:41.065 --> 00:42:42.440 align:middle line:84% It gets a little more interesting 00:42:42.440 --> 00:42:45.950 align:middle line:84% when you go to adding this next C. 00:42:45.950 --> 00:42:49.910 align:middle line:84% In order to add this C on the horizontal axis 00:42:49.910 --> 00:42:52.280 align:middle line:84% without adding anything on the vertical axis, that 00:42:52.280 --> 00:42:53.960 align:middle line:90% means that you've got an indel. 00:42:53.960 --> 00:42:57.110 align:middle line:84% And that means that you've got your A-A match. 00:42:57.110 --> 00:43:01.040 align:middle line:84% But now, to add this C, you've got a negative 2-- a penalty 00:43:01.040 --> 00:43:02.810 align:middle line:90% of minus 2. 00:43:02.810 --> 00:43:06.110 align:middle line:84% And so for the net result is a minus 1. 00:43:06.110 --> 00:43:08.990 align:middle line:84% And then, for each subsequent one, 00:43:08.990 --> 00:43:12.440 align:middle line:84% it's assumed that the extension is 00:43:12.440 --> 00:43:16.310 align:middle line:84% the same as the initial indel, which is all negative 2. 00:43:16.310 --> 00:43:19.370 align:middle line:84% And so this is an A-A match followed 00:43:19.370 --> 00:43:22.850 align:middle line:84% by one insertion, two insertions, or three 00:43:22.850 --> 00:43:25.120 align:middle line:90% insertions. 00:43:25.120 --> 00:43:27.610 align:middle line:84% And three insertions, of course, that gives you minus 5. 00:43:27.610 --> 00:43:30.790 align:middle line:84% And you just keep walking through this. 00:43:30.790 --> 00:43:35.620 align:middle line:84% Each one of these squares, essentially, 00:43:35.620 --> 00:43:39.170 align:middle line:84% is the maximum of three possibilities. 00:43:39.170 --> 00:43:43.450 align:middle line:84% The diagonal, if you follow the little yellow diagonal 00:43:43.450 --> 00:43:46.390 align:middle line:84% line from the 1 to the 0, that means 00:43:46.390 --> 00:43:50.170 align:middle line:84% you've taken an A-A match and a C-T mismatch, 00:43:50.170 --> 00:43:52.720 align:middle line:84% and the negative 1 is canceled out to positive 1, 00:43:52.720 --> 00:43:53.830 align:middle line:90% and you get a 0. 00:43:53.830 --> 00:43:56.890 align:middle line:84% Alternatively, that 0 is actually 00:43:56.890 --> 00:44:01.690 align:middle line:84% the maximum of that individual score compared 00:44:01.690 --> 00:44:11.980 align:middle line:84% to an indel from a minus 1 plus this mismatch, which is not 00:44:11.980 --> 00:44:14.650 align:middle line:84% going to be better than 0, and a minus 1 00:44:14.650 --> 00:44:18.517 align:middle line:84% plus the mismatch coming in horizontally, which is also not 00:44:18.517 --> 00:44:19.600 align:middle line:90% going to be better than 0. 00:44:19.600 --> 00:44:22.450 align:middle line:84% So you end up with 0, which is the perfect match 00:44:22.450 --> 00:44:23.710 align:middle line:90% plus a mismatch, no indels. 00:44:23.710 --> 00:44:29.080 align:middle line:90% 00:44:29.080 --> 00:44:32.910 align:middle line:84% And similarly, you can fill up the entire table this way. 00:44:32.910 --> 00:44:35.680 align:middle line:90% 00:44:35.680 --> 00:44:43.210 align:middle line:84% Finally, now, you can trace what the best scores are 00:44:43.210 --> 00:44:47.380 align:middle line:84% going from end to end here, going 00:44:47.380 --> 00:44:51.760 align:middle line:84% all the way from your A-A terminal match at the left end 00:44:51.760 --> 00:44:54.130 align:middle line:84% to the A-A terminal match at the right end. 00:44:54.130 --> 00:44:59.920 align:middle line:84% And you can see the best traceback route is going 00:44:59.920 --> 00:45:02.200 align:middle line:90% through the diagonal here. 00:45:02.200 --> 00:45:07.730 align:middle line:84% This 0 is the maximum of three possibilities-- left-right, 00:45:07.730 --> 00:45:09.700 align:middle line:90% up-down, and diagonal. 00:45:09.700 --> 00:45:14.080 align:middle line:84% Similarly, this minus 1 was the best of three possibilities. 00:45:14.080 --> 00:45:15.580 align:middle line:84% Remember, this is a global alignment 00:45:15.580 --> 00:45:17.950 align:middle line:90% that allows negative values. 00:45:17.950 --> 00:45:22.810 align:middle line:84% And its maximum was along the diagonal and so forth. 00:45:22.810 --> 00:45:26.980 align:middle line:90% 00:45:26.980 --> 00:45:31.780 align:middle line:84% That's an example of the two basic steps. 00:45:31.780 --> 00:45:34.390 align:middle line:90% You set up the scoring metrics. 00:45:34.390 --> 00:45:40.420 align:middle line:84% You set up this n by n or n by m matrix. 00:45:40.420 --> 00:45:41.770 align:middle line:90% And then you just fill it up. 00:45:41.770 --> 00:45:43.390 align:middle line:84% And that's an n squared operation. 00:45:43.390 --> 00:45:46.030 align:middle line:84% It just goes up with the number of-- 00:45:46.030 --> 00:45:47.365 align:middle line:90% the length of the two sequences. 00:45:47.365 --> 00:45:48.775 align:middle line:90% STUDENT: [INAUDIBLE]. 00:45:48.775 --> 00:45:49.775 align:middle line:90% GEORGE CHURCH: Question. 00:45:49.775 --> 00:45:52.300 align:middle line:84% STUDENT: A diagonal is not the only optimal-- 00:45:52.300 --> 00:45:53.320 align:middle line:90% GEORGE CHURCH: It's not. 00:45:53.320 --> 00:45:54.130 align:middle line:90% That's true. 00:45:54.130 --> 00:45:56.443 align:middle line:90% STUDENT: [INAUDIBLE] 00:45:56.443 --> 00:45:57.610 align:middle line:90% GEORGE CHURCH: That's right. 00:45:57.610 --> 00:46:01.270 align:middle line:84% If you have an off-diagonal that's equivalent, 00:46:01.270 --> 00:46:04.690 align:middle line:84% then that's another valid sequence alignment. 00:46:04.690 --> 00:46:07.450 align:middle line:84% And actually, it comes up quite commonly, 00:46:07.450 --> 00:46:09.580 align:middle line:84% both in global and local alignments. 00:46:09.580 --> 00:46:13.510 align:middle line:90% 00:46:13.510 --> 00:46:19.990 align:middle line:84% And then the lower right-hand corner of slide 20 00:46:19.990 --> 00:46:22.870 align:middle line:84% shows the specific interpretation 00:46:22.870 --> 00:46:24.760 align:middle line:84% of this [? Brown ?] set of errors, 00:46:24.760 --> 00:46:29.930 align:middle line:84% the particular traceback that we chose to highlight here, 00:46:29.930 --> 00:46:33.010 align:middle line:90% which is not the only one. 00:46:33.010 --> 00:46:37.990 align:middle line:84% And that's interpreted in this the same way 00:46:37.990 --> 00:46:41.600 align:middle line:90% that we interpreted the-- 00:46:41.600 --> 00:46:45.200 align:middle line:84% symbolically in an earlier slide. 00:46:45.200 --> 00:46:50.060 align:middle line:84% Now, this is also from a much earlier slide. 00:46:50.060 --> 00:46:53.240 align:middle line:84% This is the one where we had the motif 00:46:53.240 --> 00:46:55.580 align:middle line:84% to illustrate the local alignment. 00:46:55.580 --> 00:47:01.370 align:middle line:84% And on the left-hand side matrix is for local alignment using 00:47:01.370 --> 00:47:03.110 align:middle line:90% a Smith and Waterman algorithm. 00:47:03.110 --> 00:47:06.530 align:middle line:84% And on the right-hand side is the global alignment 00:47:06.530 --> 00:47:09.230 align:middle line:84% using the Needleman-Wunsch type algorithm, which we just 00:47:09.230 --> 00:47:11.040 align:middle line:90% used on a shorter sequence. 00:47:11.040 --> 00:47:22.530 align:middle line:84% And here, we've emphasized the diagonal, which 00:47:22.530 --> 00:47:25.530 align:middle line:84% gives a score of 2 and has a traceback along the magenta 00:47:25.530 --> 00:47:29.910 align:middle line:84% diagonal and would have the interpretation 00:47:29.910 --> 00:47:36.390 align:middle line:84% of the top sequence directly over the bottom sequence. 00:47:36.390 --> 00:47:38.940 align:middle line:84% On the other hand, if we look for local alignments 00:47:38.940 --> 00:47:42.120 align:middle line:84% and we do not penalize the offsets or the indels, 00:47:42.120 --> 00:47:43.550 align:middle line:90% then you can get an example. 00:47:43.550 --> 00:47:48.780 align:middle line:84% And here's another magenta traceback where we've gone-- 00:47:48.780 --> 00:47:56.240 align:middle line:84% the A-A match is not on the diagonal 00:47:56.240 --> 00:47:57.920 align:middle line:84% for the global sequence alignment, 00:47:57.920 --> 00:47:59.990 align:middle line:90% but it hasn't been penalized. 00:47:59.990 --> 00:48:03.050 align:middle line:84% So it picks up the 0's from the frame 00:48:03.050 --> 00:48:09.140 align:middle line:84% boundary cells and just picks up the positive 1 perfect match. 00:48:09.140 --> 00:48:11.330 align:middle line:84% And then, when you add a C, it picks up another one, 00:48:11.330 --> 00:48:13.010 align:middle line:90% and another C, another a. 00:48:13.010 --> 00:48:16.220 align:middle line:90% And all four add up to 4. 00:48:16.220 --> 00:48:19.910 align:middle line:84% Adding an additional base, however, does not help. 00:48:19.910 --> 00:48:24.720 align:middle line:84% Because it has to be a mismatch or an insertion or deletion. 00:48:24.720 --> 00:48:28.370 align:middle line:84% So you going from the 4 to an indel causes it to drop 00:48:28.370 --> 00:48:30.950 align:middle line:90% by minus 2, giving the two 2's. 00:48:30.950 --> 00:48:33.630 align:middle line:84% And going along the diagonal picks up a mismatch, 00:48:33.630 --> 00:48:35.330 align:middle line:90% which is a minus 1 penalty. 00:48:35.330 --> 00:48:36.630 align:middle line:90% So you just can't do better. 00:48:36.630 --> 00:48:40.130 align:middle line:84% And so that determines the edges of your local alignment. 00:48:40.130 --> 00:48:43.520 align:middle line:84% So you not only have a score and a traceback 00:48:43.520 --> 00:48:45.260 align:middle line:90% but you also have endpoints. 00:48:45.260 --> 00:48:48.250 align:middle line:90% STUDENT: [INAUDIBLE]. 00:48:48.250 --> 00:48:49.900 align:middle line:90% GEORGE CHURCH: Yes, question? 00:48:49.900 --> 00:48:50.900 align:middle line:90% You, yeah. 00:48:50.900 --> 00:48:55.180 align:middle line:84% STUDENT: Now, when this gets going, though, 00:48:55.180 --> 00:48:59.442 align:middle line:84% on the not exactly diagonal but in diagonal [INAUDIBLE]---- 00:48:59.442 --> 00:49:01.150 align:middle line:84% GEORGE CHURCH: Yes, you'll get a longer-- 00:49:01.150 --> 00:49:02.740 align:middle line:84% STUDENT: B, which is obviously not as good-- 00:49:02.740 --> 00:49:03.480 align:middle line:90% GEORGE CHURCH: And then 4. 00:49:03.480 --> 00:49:05.710 align:middle line:84% STUDENT: And then you keep going, you get back to 4. 00:49:05.710 --> 00:49:07.460 align:middle line:84% GEORGE CHURCH: And that's another example, 00:49:07.460 --> 00:49:08.585 align:middle line:90% just like the previous one. 00:49:08.585 --> 00:49:10.540 align:middle line:90% That's another valid alignment. 00:49:10.540 --> 00:49:13.510 align:middle line:90% 00:49:13.510 --> 00:49:15.850 align:middle line:90% It's still a local alignment. 00:49:15.850 --> 00:49:20.440 align:middle line:84% It doesn't have the total global endpoints. 00:49:20.440 --> 00:49:23.620 align:middle line:84% And it has an equal score to the shorter motif. 00:49:23.620 --> 00:49:25.550 align:middle line:84% STUDENT: So how do these two compare? 00:49:25.550 --> 00:49:26.420 align:middle line:90% I mean, would you-- 00:49:26.420 --> 00:49:26.930 align:middle line:90% GEORGE CHURCH: They're equal. 00:49:26.930 --> 00:49:29.430 align:middle line:84% STUDENT: --say because the other one is a longer [INAUDIBLE] 00:49:29.430 --> 00:49:32.010 align:middle line:84% GEORGE CHURCH: Well, in this case, 00:49:32.010 --> 00:49:34.940 align:middle line:84% the scoring algorithm was set up in such a way 00:49:34.940 --> 00:49:37.590 align:middle line:84% that length didn't factor into it, 00:49:37.590 --> 00:49:41.000 align:middle line:84% other than the fact that longer sequences have more chances 00:49:41.000 --> 00:49:42.260 align:middle line:90% to have more perfect matches. 00:49:42.260 --> 00:49:44.542 align:middle line:84% So in this case, they would be equivalent. 00:49:44.542 --> 00:49:47.310 align:middle line:90% 00:49:47.310 --> 00:49:54.900 align:middle line:84% As you get to more detailed substitution matrices, 00:49:54.900 --> 00:49:58.260 align:middle line:84% the chances of getting two identical scores are weaker. 00:49:58.260 --> 00:50:00.653 align:middle line:84% But with nucleic acids with this kind of simple scoring, 00:50:00.653 --> 00:50:01.695 align:middle line:90% it comes up all the time. 00:50:01.695 --> 00:50:05.010 align:middle line:90% 00:50:05.010 --> 00:50:07.320 align:middle line:90% So that's fine. 00:50:07.320 --> 00:50:14.610 align:middle line:84% We've now bargained our way down from a horrible exponential 00:50:14.610 --> 00:50:17.280 align:middle line:84% potential way of doing alignments 00:50:17.280 --> 00:50:21.120 align:middle line:84% to something which scales by n times 00:50:21.120 --> 00:50:26.910 align:middle line:84% m where those are the lengths of the two sequences. 00:50:26.910 --> 00:50:30.660 align:middle line:84% You could think of this as a rectangular matrix 00:50:30.660 --> 00:50:32.640 align:middle line:84% such as the ones we've been doing. 00:50:32.640 --> 00:50:38.280 align:middle line:84% And both the time and the space or memory requirements 00:50:38.280 --> 00:50:45.030 align:middle line:84% for the algorithm will scale by this quadratic relationship. 00:50:45.030 --> 00:50:49.640 align:middle line:84% And the amount of time in memory is modest. 00:50:49.640 --> 00:50:52.330 align:middle line:90% 00:50:52.330 --> 00:50:55.680 align:middle line:84% So in absolute terms, it would be on the order of one 00:50:55.680 --> 00:50:58.920 align:middle line:84% comparison-- that's that maximum comparison-- 00:50:58.920 --> 00:51:02.130 align:middle line:90% 00:51:02.130 --> 00:51:05.530 align:middle line:84% three multiplication steps and in computing the entry. 00:51:05.530 --> 00:51:09.420 align:middle line:84% And the memory could be on the order of a byte. 00:51:09.420 --> 00:51:12.940 align:middle line:84% Data structure could be integer, or it could be floating point. 00:51:12.940 --> 00:51:16.830 align:middle line:84% And again, you have to have some way of finding 00:51:16.830 --> 00:51:19.685 align:middle line:90% the entries in the table. 00:51:19.685 --> 00:51:20.310 align:middle line:90% So that's fine. 00:51:20.310 --> 00:51:21.570 align:middle line:90% It scales gracefully. 00:51:21.570 --> 00:51:23.730 align:middle line:90% But how big is it? 00:51:23.730 --> 00:51:28.440 align:middle line:84% Let's say we had two megabase genomes. 00:51:28.440 --> 00:51:31.890 align:middle line:84% In order to find entries of that size, 00:51:31.890 --> 00:51:34.290 align:middle line:84% you might want to set aside 4 bytes. 00:51:34.290 --> 00:51:40.320 align:middle line:84% And so you have the 4 bytes times 10 to the sixth squared. 00:51:40.320 --> 00:51:42.083 align:middle line:90% Or this is just ballpark. 00:51:42.083 --> 00:51:44.250 align:middle line:84% There are various ways you could squeeze this a bit. 00:51:44.250 --> 00:51:47.280 align:middle line:84% But this is 4 gigabytes of memory. 00:51:47.280 --> 00:51:53.280 align:middle line:84% And for a gigahertz CPU, you might 00:51:53.280 --> 00:51:58.370 align:middle line:84% be able to do a million entries per second so that with 10 00:51:58.370 --> 00:52:01.400 align:middle line:84% to the 6 squared entries, that's about 10 days. 00:52:01.400 --> 00:52:03.860 align:middle line:84% Now, that's a fairly small genome. 00:52:03.860 --> 00:52:05.990 align:middle line:84% Most genomes are bigger than 1 megabase. 00:52:05.990 --> 00:52:07.730 align:middle line:84% And so when we had the discussions 00:52:07.730 --> 00:52:10.310 align:middle line:84% at the beginning of the Genome Project, 00:52:10.310 --> 00:52:12.650 align:middle line:84% one of the things the computer people brought up 00:52:12.650 --> 00:52:19.130 align:middle line:84% was, how are we going to compare a billion base pairs 00:52:19.130 --> 00:52:21.590 align:middle line:84% with a billion base pairs if the goal of this project 00:52:21.590 --> 00:52:24.650 align:middle line:84% is to do the three billion base pair human genome? 00:52:24.650 --> 00:52:30.770 align:middle line:84% And of course, back then, most computers were 4 gigabytes, 00:52:30.770 --> 00:52:37.340 align:middle line:84% and a gigahertz was a quite remarkable computer. 00:52:37.340 --> 00:52:39.320 align:middle line:84% And of course, the answer was that we 00:52:39.320 --> 00:52:42.050 align:middle line:84% weren't going to do a full dynamic programming 00:52:42.050 --> 00:52:44.300 align:middle line:84% of the human genome against itself. 00:52:44.300 --> 00:52:48.200 align:middle line:84% We were going to cheat in various ways. 00:52:48.200 --> 00:52:55.200 align:middle line:84% And it took the recognition that it really 00:52:55.200 --> 00:52:59.190 align:middle line:84% was practical and, biologically, not much of a shortcut 00:52:59.190 --> 00:53:03.840 align:middle line:84% to look for little anchor points that would tell you 00:53:03.840 --> 00:53:06.970 align:middle line:84% that maybe the sequences don't align end to end, 00:53:06.970 --> 00:53:08.820 align:middle line:84% but there's some anchor point where you have 00:53:08.820 --> 00:53:11.640 align:middle line:84% enough bases or enough amino acids in a row 00:53:11.640 --> 00:53:13.680 align:middle line:84% that you can say, OK, here's one point where 00:53:13.680 --> 00:53:14.970 align:middle line:90% they definitely align. 00:53:14.970 --> 00:53:17.370 align:middle line:84% Let's now make reasonable assumptions 00:53:17.370 --> 00:53:23.670 align:middle line:84% by how many indels there can be, for example, 00:53:23.670 --> 00:53:27.330 align:middle line:84% by knowing how different the two sequences are. 00:53:27.330 --> 00:53:30.420 align:middle line:84% And so if you know the differences of sequences, 00:53:30.420 --> 00:53:33.840 align:middle line:84% then you can say, I'm not going to allow 00:53:33.840 --> 00:53:36.210 align:middle line:84% more than a reasonable number of indels 00:53:36.210 --> 00:53:38.130 align:middle line:84% based on how different the sequences are. 00:53:38.130 --> 00:53:40.920 align:middle line:84% And you make a band which is a narrow width-- 00:53:40.920 --> 00:53:44.220 align:middle line:84% here's a fairly extreme example where we have a width of 3. 00:53:44.220 --> 00:53:46.080 align:middle line:84% And so rather than doing a full n 00:53:46.080 --> 00:53:50.820 align:middle line:84% squared matrix where you filled up the entire thing, 00:53:50.820 --> 00:53:55.350 align:middle line:84% we just do this band, which is on the order 00:53:55.350 --> 00:54:00.090 align:middle line:84% of the width of the band times the length of the longest 00:54:00.090 --> 00:54:02.220 align:middle line:90% sequence. 00:54:02.220 --> 00:54:06.330 align:middle line:84% Now, this doesn't look very impressive for this case 00:54:06.330 --> 00:54:09.554 align:middle line:84% because n is small, and w is relatively large. 00:54:09.554 --> 00:54:14.220 align:middle line:84% But if n were billions and w were, say, 3 to 5, 00:54:14.220 --> 00:54:16.990 align:middle line:84% then it would be a very significant savings. 00:54:16.990 --> 00:54:18.550 align:middle line:90% So there's two key things here. 00:54:18.550 --> 00:54:20.580 align:middle line:84% One is the banding, and the other 00:54:20.580 --> 00:54:23.110 align:middle line:90% is getting the anchor points. 00:54:23.110 --> 00:54:28.080 align:middle line:84% So summary for this half of the talk 00:54:28.080 --> 00:54:33.810 align:middle line:84% is that dynamic programming is really the rigorous way 00:54:33.810 --> 00:54:35.055 align:middle line:90% to compare two sequences. 00:54:35.055 --> 00:54:36.840 align:middle line:84% And after the break, we'll see how you 00:54:36.840 --> 00:54:39.210 align:middle line:90% can compare multiple sequences. 00:54:39.210 --> 00:54:42.240 align:middle line:84% We need to work towards a statistical interpretation 00:54:42.240 --> 00:54:43.260 align:middle line:90% of these alignments. 00:54:43.260 --> 00:54:47.490 align:middle line:84% That's going to require some test sets-- 00:54:47.490 --> 00:54:49.150 align:middle line:90% sorry, some training sets-- 00:54:49.150 --> 00:54:51.030 align:middle line:84% where you can see how it actually 00:54:51.030 --> 00:54:57.240 align:middle line:84% behaves on real biological populations of sequence 00:54:57.240 --> 00:54:58.620 align:middle line:90% alignments. 00:54:58.620 --> 00:55:02.357 align:middle line:84% We need to compute either a global or a local alignment. 00:55:02.357 --> 00:55:04.440 align:middle line:84% And you've seen algorithms for doing each of those 00:55:04.440 --> 00:55:07.260 align:middle line:84% and how there's important but subtle differences 00:55:07.260 --> 00:55:08.322 align:middle line:90% between them. 00:55:08.322 --> 00:55:09.780 align:middle line:84% And we've talked about ways that we 00:55:09.780 --> 00:55:17.700 align:middle line:84% can improve the algorithm tremendously using 00:55:17.700 --> 00:55:21.075 align:middle line:84% the simple scoring functions or complicated ones that 00:55:21.075 --> 00:55:22.200 align:middle line:90% are determined empirically. 00:55:22.200 --> 00:55:23.580 align:middle line:90% So let's take a little break. 00:55:23.580 --> 00:55:27.500 align:middle line:84% And we'll come back and talk about multisequence alignment. 00:55:27.500 --> 00:55:28.000 align:middle line:90%