Practical: Matching real reads

Loading...

Do curso por Johns Hopkins University

Algoritmos para Sequenciamento de DNA

367 classificações

Johns Hopkins University

Algoritmos para Sequenciamento de DNA

367 classificações

Curso 4 de 8 em Specialization Genomic Data Science

We will learn computational methods -- algorithms and data structures -- for analyzing DNA sequencing data. We will learn a little about DNA, genomics, and how DNA sequencing is used. We will use Python to implement key algorithms and data structures and to analyze real genomes and DNA sequencing datasets.

Na lição

DNA sequencing, strings and matching

This module we begin our exploration of algorithms for analyzing DNA sequencing data. We'll discuss DNA sequencing technology, its past and present, and how it works.

Module 1 Introduction1:32

Lecture: Why study this?4:29

Lecture: DNA sequencing past and present3:08

Lecture: Genomes as strings, reads as substrings5:02

Lecture: String definitions and Python examples3:33

Practical: String basics 7:14

Practical: Manipulating DNA strings 7:00

Practical: Downloading and parsing a genome 6:29

Lecture: How DNA gets copied3:37

Optional lecture: How second-generation sequencers work 7:10

Optional lecture: Sequencing errors and base qualities 6:37

Lecture: Sequencing reads in FASTQ format4:39

Practical: Working with sequencing reads 11:11

Practical: Analyzing reads by position 6:58

Lecture: Sequencers give pieces to genomic puzzles5:41

Lecture: Read alignment and why it's hard3:09

Lecture: Naive exact matching10:18

Practical: Matching artificial reads 6:53

Practical: Matching real reads 7:24

Conheça os instrutores

Ben Langmead, PhD
Assistant Professor
Computer Science
Jacob Pritt

Department of Computer Science

Now that we've implements native exact matching and use it to match artificial

reads against the genome, let's do the same thing with real sequencing reads.

So we've already downloaded the 5X genome.

The other link we've given you is to a set of actual

sequencing reads from that same genome.

So you can download that, and I'm going to do that with Wget.

And to read the FASTQ,

I'm going to use the same readFastq function that we developed in Practical 4.

So I'm just going to paste that in here.

And I'm going to I'm going to use that to actually read the FASTQ file.

And this function actually returns qualities too,

which in this case we're not interested in.

So, I'm just going to use a blank to indicate that we don't want to save that

to a variable.

I'm going to paste in the name of this file.

As this has downloaded the sequencing reads from that file and

that's a set of a 1000 reads.

So, now let's do the same thing we did above with artificial reads

in practical six and actually match these, align these against the genome and

see how many of them match.

So, I'm going to keep track of how many reads have matched.

So for each read r, in my list of reads.

I'm going to get the set of matches from naive.

And n is just a count of the total number of reads that we've processed.

So now I can say, if the length of the set of matches returned is greater than zero,

that means that they read aligned in at least one place, so we found a match.

So I can increment numMatched.

And now I'm just going to print the proportion of reads that matched

correctly.

So in order to run that, starting the k status running, so

it's taking a little while.

And only seven out of 1,000 reads match the genome.

>> Okay, so that might be surprising.

These are sequencing reads actually from this genome, and yet

only 7 out of 1,000 of them actually match the genome.

So why would that be?

Well, we actually learned a couple of reasons why.

A sequence obtained with a DNA sequencer might not match

back to the genome exactly.

One reason is because of sequencing errors, so

it could be that some of the bases were read incorrectly by the sequencer.

And then another reason is because there might be differences between

the individual that was sequenced,

in this case it's a virus, and the reference genome sequence.

So in this case it's probably mostly the sequencing errors that are the reason why

many of these reads are not aligning.

So one of the things that we're going to try next is to just take

part of the read, right?

It could be that the whole read, all 100 bases, don't necessarily match.

But if we just take a piece, like maybe the first 30 bases or so.

That has a greater chance of matching exactly.

All right, so we'll try that next.

>> Okay. So, I just can make a little change to get

the first 30 bases.

We'll set r equal to, just the 30 base prefix of r.

I'm going to run it again.

And this time, we're 459 out of 1000 reads match.

So, a little under one half.

>> So this still isn't very close to a thousand.

This is less than half of the reads aligning to our reference genome.

So you might be able to think of another reason why maybe not so

many of the reads are aligning.

In this case its because the genome is double stranded and

the reads can come from one strand or they can come from the other.

But the exact matching problem we've set up is only looking on

one of those two strands.

It's not looking on both.

So that could explain in large part why we're seeing

less than half of the reads actually matching back to the genome.

So the next thing to try would be not just to try to match the read to the genome,

but also to try to match the reverse complement of the read to the genome.

That will cover the case where the read comes from the other strand of the genome.

>> Yeah, so I can do that using the reverse compliment function,

which we developed in Practical 2.

I'm just going to paste this in here.

Now I'm going to use basically the same code from up here.

But I'm going to do two steps here.

So I'm setting matches to be the naive exact match between r and the genome.

But I'm going to do something similar, I'm going to say matches.extend

with naive of the reverse compliment

of r in the genome.

So what dot extend means is it just takes that list and

adds on to it the results of another list.

So the first line here, matches equals naive of r in the genome

is going to test whether that read matches in the forward direction on the genome.

In the next line, we'll add onto it the results

of any matches of that read to the reverse complement of the genome.

So I'm going to run this now.

>> Oh [LAUGH].

So the problem here is that our.

Our dictionary in the reverse complement function doesn't allow for

there to potentially be an n in a sequencing read.

>> Right. >> And yet, like we said before,

sometimes a sequencing read has an n meaning no competence.

So we just have to add one more key value pair to our dictionary there.

>> So I'll just add the base n and this will just match to an n.

Since in both cases we don't know what it is.

>> Yeah. >> So I mean update that and run it again.

>> It's taking a little longer because it has to do more matching now.

>> There we go 932 out 1,000 reads have matched.

>> Yeah. So still not 100% But

this is to be expected, and because of sequencing errors and things like that.

But once we figured out that we had missed one whole strand of the genome, and

once we took a prefix of the read instead of the entire read,

we were able to match nearly all the reads in our data set.

But the fact that we had to take the first 30 base pairs or so

is a pretty good indication that exact matching

is not necessarily going to be what we want to do, right?

We want to allow for sequencing errors,

we want to allow for differences between individuals.

So that's why later in the course we'll

examine some approximate matching algorithms.

>> Yeah.

Explore nosso catálogo

Registre-se gratuitamente e obtenha recomendações, atualizações e ofertas personalizadas.

O Coursera proporciona acesso universal à melhor educação do mundo fazendo parcerias com as melhores universidades e organizações para oferecer cursos on-line.

© 2019 Coursera Inc. Todos os direitos reservados.

Baixar na App Store

Baixar no Google Play