SEQanswers

Go Back   SEQanswers > Applications Forums > De novo discovery

Similar Threads
Thread Thread Starter Forum Replies Last Post
RNA-seq SNP-calling without a complete reference shoegame2001 RNA Sequencing 6 07-04-2012 01:55 AM
How to generate de novo assembly sequence from complete genomics data? ssnowfox Bioinformatics 2 04-19-2012 10:34 PM
SNP calling from a reference sequence blackrabite Genomic Resequencing 2 05-21-2011 09:48 PM
Hierarchical reference-free SNP calling Marius Bioinformatics 1 12-27-2010 09:38 AM
PubMed: Direct resequencing of the complete ERBB2 coding sequence reveals an absence Newsbot! Literature Watch 0 04-18-2008 07:10 AM

Reply
 
Thread Tools
Old 02-16-2012, 04:40 AM   #1
fcr
Member
 
Location: Seville

Join Date: Jan 2012
Posts: 19
Default De novo SNP calling in absence of complete reference assembly

Dear all,

I would like to call SNPs on a diploid genome in an absence of a reference genome or complete assembly. I know that is possible to do this with cortex but has a lower sensitivity than mapping-based approaches. (And also requires around 30 X coverage).

The ultimate assembly for our species is not ready yet. We do have access to about 1 million scaffolds at this moment. Additional WGS reads with 20x coverage are available for several individuals. I plan to start the SNP calling now. My idea is to map the reads against the scaffolds and then use FreeBayes for calling. what do you think? Would 1 SAMtools work with 1 million scaffolds?

I would really appreciate receiving any advice or comments.

Cheers,
fcr

Last edited by fcr; 02-16-2012 at 04:49 AM. Reason: incomplete post
fcr is offline   Reply With Quote
Old 02-16-2012, 05:04 AM   #2
lh3
Senior Member
 
Location: Boston

Join Date: Feb 2008
Posts: 693
Default

I think you should map your reads to the assembly and then do SNP calling. SAMtools should in principle work, but I have not tried.
lh3 is offline   Reply With Quote
Old 02-16-2012, 05:17 AM   #3
Zam
Member
 
Location: Oxford

Join Date: Apr 2010
Posts: 51
Default

Re Cortex:
1. You have much more than 30x coverage if you have many samples at 20x
2. It's not as simple as "you need 30x" for Cortex. But you are absolutely right that an assembly approach will be less sensitive to SNPs.


Re what to do
- it depends what you want to achieve. Do you want a conservative small set of SNPs for building a genetic map, or a big sensitive set for some other purpose etc.

If you have the time, then try both methods (mapping/assembly) and compare. If you are doing population genetic studies, then experience suggests that you will need to be very careful with SNP calls based on an assembly that is not high quality, as it is easy for assembly artefacts to look like interesting scientific finds in your SNPs.
Zam is offline   Reply With Quote
Old 02-16-2012, 07:10 AM   #4
fcr
Member
 
Location: Seville

Join Date: Jan 2012
Posts: 19
Default

Hi Zam,

Thanks a lot.

Cortex:
1. True, the distribution of coverage will include regions above 30x.
2. What are the Computational needs for 10 individuals with 2.9 Gbp genome? In your server you stated "10 humans on a 256Gb RAM server" How long this takes? Would it be possible to call SNPs with less RAM?

What to do:
This is a 60 X coverage genome. I would assume that many of the scaffolds are bona fide and that many of the changes (adding more libraries) are going to affect mainly to the connection among scaffolds rather than disrupting them...but I might be wrong and shouldn't guess. The main interest are; 1. develop genome-wide set of markers and 2. do some population inferences by estimating Fst, Pi and Ne.

So you think is too risky using scaffolds?


Cheers,
Fernando
fcr is offline   Reply With Quote
Old 02-16-2012, 07:26 AM   #5
Zam
Member
 
Location: Oxford

Join Date: Apr 2010
Posts: 51
Default

Hi Fernando


>True, the distribution of coverage will include regions above 30x.
One of the examples in our paper is of SNP calling in 10 samples each sampled to 6x,
for example.

2. Actually, you could call on 10 individuals with much less than 256Gb of RAM. You need 256Gb to hold all of ALL of their genomes at the same time. But lots of the genome is either monomorphic, or doesn't consist of things Cortex can call. So you could do those 10 samples in ~80Gb of RAM (for comparison I've just done 85 humans in 320 Gb of RAM).
The trick is to call on the joint graph (1 colour, probably needs 80Gb RAM) and then pull out just the variants and make a graph just of the variants. Then "multicolourise" the graph and make a 10-colour graph of the variants only, and genotype everyone in that.
Uses far less memory.

How much coverage do your 10 samples have? Is the 60x individual a different sample?

I'm not saying it is too risky with scaffolds, just that if you find something really exciting, you need to do some work making sure it's not an artefact. I've seen people have to work very hard to avoid problems with the chimp genome.

best

Zam
Zam is offline   Reply With Quote
Old 02-16-2012, 07:33 AM   #6
lh3
Senior Member
 
Location: Boston

Join Date: Feb 2008
Posts: 693
Default

With 60X, you should be able to get an assembly decent enough for most analyses. This is true for human. Nonetheless, Zam is right that misassembly may cause artifacts. You have to live with it. If you are careful enough, you can greatly reduce the effect of that. Also beware that there will be reference bias when estimating population statistics (i.e. individuals closer to the reference will be mapped better).
lh3 is offline   Reply With Quote
Old 02-16-2012, 07:37 AM   #7
Zam
Member
 
Location: Oxford

Join Date: Apr 2010
Posts: 51
Default

Just to clarify one thing (and agree with Heng) - my understanding is that Fernando doesnt want to have to wait until his assembly is finished (I mean done/completed, not finished by manual finishers), and wants to get on with it and start calling now. That's what got me nervous about artefacts.
Zam is offline   Reply With Quote
Old 02-17-2012, 08:56 AM   #8
fcr
Member
 
Location: Seville

Join Date: Jan 2012
Posts: 19
Default

Hi,

Yes, Zam got it right. I want to start calling SNPs now. The assembly is unfinished and it's going to take time polishing it (~1000,000 scaffolds now). In response to Zam, the assembly is based on an individual, and the estimated coverage is 60X.

The other 10 individuals have 20 X coverage and i want to use them for SNP calling and perhaps "pilot" genotype calling. I think is worthy advance on this, even if in the future a second calling based on the assembly will help to verify/reject candidate regions of interest.

lh3: Thanks for your comment about the reference bias when estimating the population statistics...I will keep that in mind.

Cheers,
Fernando
fcr is offline   Reply With Quote
Old 09-21-2012, 02:35 AM   #9
rururara
Member
 
Location: montreal

Join Date: Jan 2011
Posts: 31
Default De novo SNP calling in absence of complete reference assembly

Hai all,

What about if the incomplete reference genome like papaya? The available information on papaya are scaffolds and contigs. Is it possible to use papaya scaffolds as a reference to align against my reads? In my case, the objective is to discover the SNPs.
rururara is offline   Reply With Quote
Old 09-21-2012, 02:49 AM   #10
Zam
Member
 
Location: Oxford

Join Date: Apr 2010
Posts: 51
Default

Hi Rururara
Are you working on the same project as Fernando or a different one? If different, how many samples are you trying to discover SNPs in, and what are their depths od coverage and with what technology. Finally, sorry for ignorance, but what is the ploidy of papaya?
regards
Zam
Zam is offline   Reply With Quote
Old 09-21-2012, 02:57 AM   #11
fcr
Member
 
Location: Seville

Join Date: Jan 2012
Posts: 19
Default

Hi Zam,

Rururara is not working in the same project as me. If papaya is a diploid, he could probably use the papaya scaffolds with the "Coordinates Only" option during the calling with cortex_var (actually a acompanying script called runcalls.pl). Right?

Cheers,
Fernando
fcr is offline   Reply With Quote
Old 09-21-2012, 03:06 AM   #12
Zam
Member
 
Location: Oxford

Join Date: Apr 2010
Posts: 51
Default

Yes, and to explain that in more detail:
Rururura:

1.If you have one diploid sample you can de novo discover variants using Cortex, and then use your contigs/scaffolds to assign them coordinates. This is what Fernando meant by "CoordinatesOnly", an option for Cortex's new wrapper script.

2. If you have several samples, then you can do two things
a) You can also use the Cortex "population filter" to classify putative variants as repeat/error/polymorphism - this method is robust to reference assembly errors - it catching missing collapsed repeats in the reference - and this will give you a high quality set of variants
b) you could use this method to look into the quality of the reference and annotate regions which you trust and do not trust.

Zam
Zam is offline   Reply With Quote
Old 09-21-2012, 03:11 AM   #13
rururara
Member
 
Location: montreal

Join Date: Jan 2011
Posts: 31
Default

Hi Zam & fcr,

Yup, we are not in the same team. Hehe. Papaya is diploid. I have 3 samples and one of the sample is parental lines. I'm not sure yet the depth coverage as I am still not getting any sequencing information from the company, but soon I will. Papaya is sequence using HiSeq platform.
rururara is offline   Reply With Quote
Old 09-21-2012, 03:26 AM   #14
Zam
Member
 
Location: Oxford

Join Date: Apr 2010
Posts: 51
Default

Hi there- when you say one of the samples is parental, does that mean you have two parents and 2 F1 samples, and you have sequenced one parent and both progeny?
Zam
Zam is offline   Reply With Quote
Old 09-21-2012, 03:31 AM   #15
rururara
Member
 
Location: montreal

Join Date: Jan 2011
Posts: 31
Default

Definitely yes. Is there any concern about that? Do u mind to share? Anyway, I would like to try this approach whereby I assemble the parental reads with scaffold and use it as a reference sequence to align against the other two progeny. What do u think?
rururara is offline   Reply With Quote
Old 09-21-2012, 03:34 AM   #16
Zam
Member
 
Location: Oxford

Join Date: Apr 2010
Posts: 51
Default

It's fine. I have more questions - shall we take it offline? I'm zam AT well.ox.ac.uk
Zam is offline   Reply With Quote
Reply

Thread Tools

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off




All times are GMT -8. The time now is 07:02 PM.


Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2021, vBulletin Solutions, Inc.
Single Sign On provided by vBSSO