Jellyfish kmer counting diploid

[Bourney]

Hi there,

I've used jellyfish kmer count to estimate the genome size of a diploid organism. I get double peaks, which are expected, and the first one is twice the size, which is also expected when sequencing diploid. However, I've always been told the protocol is to disregard the first one, as it represents the heterozygous alleles, and use the peak of the second one for genome size estimation. But this advice has mainly been for haploid sequencing.
And this paper, https://academic.oup.com/dnaresearch...inbred-line-of, looking at a tetraploid, used the first peak in their estimation.

Would I be justified in using the first peak for estimation of a diploid genome? I can't think why I wouldn't, but don't want to have overlooked something critical.

For what it's worth, when using the first peak, I get an estimation of 170 Mbp, which is around the right size (same as sister-species), and when using SOAPdenovo and Abyss to assemble, I get an assembly twice that at 330-340 Mbp. I'm assuming that the assemblers aren't handling the diploidy overly well. So, as an aside, does anyone have any advice on making the assemblers incorporate the heterozygosity of a diploid genome?

Cheers

[Brian Bushnell]

The correct way to calculate the (haploid) genome size of a diploid is:

0.5*(# kmers in first [heterozygous] peak)+
1*(#kmers in second [homozygous] peak)+
2*(#kmers in third [2-copy] peak)+
...etc

Normally you can ignore the high peaks as they become increasingly hard to distinguish don't contribute a lot. But you absolutely can't disregard the first peak; in highly heterozygous organisms, that can constitute the bulk of the genome!

BBMap's KmerCountExact will generate a peaks file that does this automatically, when the peaks are clear (as in your example):

kmercountexact.sh in=reads.fq khist=khist.txt peaks=peaks.txt ploidy=2

A header line in peaks.txt will indicate genome size (and calculated ploidy, if you don't specify it).