Suggestion for mitochondrial genome assembly

[bioman1]

1.I have denovo assembled my illumina paired-end reads (2x101 bp with abyss genome assembler and got contigs.
2.I have blasted my denovo assembled contig assembly against my closest reference mitochondria genome (from NCBI organelle genome)
3.Extracted contigs which have greater than 80% identity with my reference mitochondria genome
4.Through blast alignment, I extracted 40 putative mitochondrial contigs from my denovo assembled genome


Now I want my extracted contigs as single mitochondrial genome assembly as in NCBI organelle genome (for e.g.NC_012119.1: http://www.ncbi.nlm.nih.gov/nuccore/...9?report=fasta)


I tried reference mapping with my close reference, iterative mapping using genious, codon code assembler and also MITObim. But could not able to get
as single complete genome. Please let me know how to get single complete genome?

[juan.isaza]

you can reoder your denovo assembled contigs into a single molecule (with some gaps) using ABACAS (it needs MUMmer) and the closest reference mitochondria genome, after that you can use IMAGE to fill some gaps with your illumina reads.

ABACAS ----> http://abacas.sourceforge.net/
IMAGE ----> https://www.sanger.ac.uk/resources/s...e/pagit/#IMAGE

[dsenalik]

Since you used ABySS, you can use ABySS-Explorer to investigate connections between contigs. It will probably be difficult, since most mitochondrial genomes contain at least some repeated sequences, and sequence shared with plastid genome, and nuclear genome, which may make the assembly ambiguous. And even if your reference is the same species, there could be major rearrangements. If it is a different species, there most certainly are rearrangements. Good luck!

[Brian Bushnell]

Your coverage may be too high for good assembly. Approximately how much coverage do you have? Once you have isolated the mitochondrial contigs, you can map to them and collect all the reads that hit. Then subsample or normalize those reads to an appropriate depth (maybe 100x) and assemble that. Error-correction may help as well.

[Matt Kearse]

Brian's suggestion sounds great although I'd tweak it to map all your reads to a combination of your contigs and any related reference mitochondria genomes you can find. And when collecting the reads that hit, make sure you include the unmapped mates of reads where only one from the pair hit. Then de novo assemble that data set, subsampling at an appropriate depth as suggested by Brian.

An alternative is that once you know what an appropriate subsample fraction for mitochondria is from the above step, de novo assemble the same fraction of your entire original data set.

[bioman1]

Thanks all for your suggestion. When I map my reads to denovo assembled genome (nuclear+mito) i get around 5x coverage.

Also, I tried by mapping my reads to reference mitochondria genome and I got mapped reads. Then mapped reads where assembled with sub-sampled WGS illumina data (10% random reduction) with MITObim. However, I got lot of gaps in the assembly.

So, I tried second method, by denovo assembling the WGS illumia paired-end reads and collected putative contigs by blasting against my closest reference mitochondiral genome.

1.I get around 40 contigs which got hit above 80% similarity and are above 1Kbp.
2.I get around 3 contigs which are above 3Kbp and 1 contig around 5kbp. Does that mean this one contig 5Kbp will be mitochondrial contig, can I select this for iterative mapping and assembling?

Also I would like to try Brian's suggestion, to sub-sample with apporpriate depth (100x), how to get the reads which have greater depth (I read some of papers that they reads which have greater depth and denovo assemble). Does your BBMap can get reads for appropriate depth? If not, can you suggest some tools?

[Brian Bushnell]

The BBTools package comes with BBNorm, which can bin reads by depth. Sometimes we use that to separate mitochondrial reads from nuclear reads. e.g:

bbnorm.sh in=reads.fq passes=1 keepall lowbindepth=20 highbindepth=40 outlow=nuclear.fq outmid=unknown.fq outhigh=mito.fq

But for this to work you need to know the coverage specifically over the mitochondrial part to figure out what the thresholds should be. You can determine that with my pileup program if you have a sam file (or generate it directly from BBMap):

pileup.sh in=mapped.sam covstats=stats.txt

That will tell you the average coverage for every contig, so if you know which are the mito contigs and which are the genomic contigs, it will tell you what limits to use.