Extending contigs N50 (SOAPdenovo)?

[condomitti]

Hellow fellows,

I've been assembling a genome of a snake with SOAPdenovo using big set of paired-end reads that were sequenced in Illumina hiScan.
I already made a lot of assemblies tuning the parameters and using filtered/non-filtered reads, but I'm always getting a low contig N50 of 1.4k (against >30k on references).

I've read something about breaking the scaffolds in contigs and using the reads again to map one end to the gap but I couldn't find a detailed explanation on how to do that. The fact is that in all assemblies that took this approach they could extend contigs lengths and N50 as a consequence.

Could someone recommend me a text to read about that? has any of you ever done this?


Also, another doubt I have is regarding evaluating scaffolds with QUAST (http://quast.bioinf.spbau.ru/manual.html)... I see it analyzes the scaffold file in two ways, and one of them gives me a bigger N50 (i.e. 3.5k) but I don't know the differences between the two results quast gives me.
No problem when evaluating contigs though.
Any clues here too?


Thanks a lot in advance!

Condomitti.

[condomitti]

No one has tips on this??

[AdrianP]

There is no different ways of calculating N50, only one way. Find a perl script that does it.

Doesn't soap already provide you with scaffold and contig files? There should be no need to break scaffolds into contigs.

Also, does SOAP not give you the stats? the n50, the longest contig...?

What you are suggesting to do with reads is called read walking, where you maps reads to the ends of contigs to extend. You will generate missassemblies due to repeats while doing that.

The best thing you can try is other assemblers, other parameters.

[condomitti]

Thank you for your reply, AdrianP!

Yes, SOAP gives me all the assembly stats. I have N50 value indeed. My concern is regarding the low N50 value.

I have tried SGA, which only increased that value in a few units.



cheers,
Condomitti.

[AdrianP]

You have illumina data, SGA might not be your best choice. A few questions.

What is the rough nucleotide coverage of the genome?
What about the genome size?
How long are your reads? And how long is the insert size?

[condomitti]

The genome size is 2.2Gbp

Considering contigs > 800bp, the nucleotide coverage is 524.993.890bp


Reads vary from 40-100 and insert size 300bp.


N50 1549bp
largest contig: 28.559bp


Cheers,
Condomitti.

[AdrianP]

The genome size is 2.2Gbp

Considering contigs > 800bp, the nucleotide coverage is 524.993.890bp


Reads vary from 40-100 and insert size 300bp.


N50 1549bp
largest contig: 28.559bp


Cheers,
Condomitti.

I think you gave me the assembly size rather than the nucleotide coverage. Nucleotide coverage is how many reads overlap any given DNA sequence from the genome.

Your reads vary in length. That is not normal for illumina sequencing, did you trim them, or do you have different libraries?

What kmer values did you use when assembling with SOAP?

[condomitti]

Ohh you are right, sorry about that...

I'm working with ~130x fold.

I did trim them, and applied some filters to remove duplications etc.

I have tried some different values for kmer, using both single and multi-kmer strategies.

With single kmer, the one that generated better results was 65.
Using multi-kmer 61-71 I could get that result I've written above.

[AdrianP]

Using untrimmed libraries, try SPAdes wither kmers 23,33,43,53,63,73

[condomitti]

Thanks AdrianP! I'll take a look.

Cheers,
Condomitti.