SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
Genome size estimation using BBMAP syintel87 Bioinformatics 12 01-04-2018 10:27 AM
Size issue of Pacbio Reads rusty1947 Pacific Biosciences 11 01-11-2015 11:59 AM
Genome size estimation for paired end data mandar.bobade60 General 4 11-19-2014 07:56 PM
Genome size estimation-jellyfish bioman1 Bioinformatics 3 08-18-2014 11:14 AM
Genome size estimation moinul De novo discovery 9 04-04-2014 03:22 AM

Reply
 
Thread Tools
Old 07-14-2016, 11:35 AM   #1
jpummil
Member
 
Location: Fayetteville, AR

Join Date: Apr 2014
Posts: 85
Default Genome Size Estimation from PacBio Raw Reads

So, working on a de novo assembly using Canu, and it seems to be VERY sensitive to the genomeSize=XXX parameter which is required. As it is a new project, no one has an actual "size" on it (checked T Ryan Gregory's site...nothing similar there either).

So, I am using BBMap suite, specifically...the "kmercountexact.sh" component. Waiting on a compute node right now with >64GB of ram to run, but have it set as follows: kmercountexact.sh in=filtered_subreads.fastq khist=khist.txt peaks=peaks.txt out=genomesize.txt

As Brian Bushnell is active on here, I was hoping to inquire about using this on PacBio specifically...anything I need to be more specific about on the options? Also, can I specify both of my PacBio files as arguments? I have both a .fastq of the long reads as well as a .fasta of much shorter reads supplied by the sequencer people. I know it can do PE files as in= and in2=, but what about to essentially "single" reads?
jpummil is offline   Reply With Quote
Old 07-14-2016, 02:14 PM   #2
Brian Bushnell
Super Moderator
 
Location: Walnut Creek, CA

Join Date: Jan 2014
Posts: 2,707
Default

Hi Jeff,

Unfortunately, I don't have a good method for this. I've tried kmercountexact, and it does not work on raw PacBio reads due to the high error rate. I do not know of a better method for genome size estimation than assembling, with Falcon, for example. Sorry!

If you have multiple files, though, you can enter them comma-delimited, like this:

Code:
kmercountexact.sh in=filtered1.fq,filtered2.fq
Not all tools support that, but Tadpole, KmerCountExact, and Dedupe do.
Brian Bushnell is offline   Reply With Quote
Old 07-14-2016, 06:48 PM   #3
jpummil
Member
 
Location: Fayetteville, AR

Join Date: Apr 2014
Posts: 85
Default

Thanks for the quick response, Brian!

Good to know about the comma-delimited method for multiple entries. Unfortunate to hear about the PacBio error issue when trying to determine genome size. I thought about this a bit and am wondering if the pre-processing Canu does to the data could be used prior to trying kmercountexact? It outputs a couple of files during its run which trim, then correct the reads:

<filename>.trimmedReads.fasta.gz

then

<filename>.correctedReads.fasta.gz

Of course, they have been processed WITH the genomeSize estimate provided at run time and I'm not certain of how much that might have influenced any trimming or correction. I might try and contact Phillippy or Koren and inquire further ;-)
jpummil is offline   Reply With Quote
Old 07-15-2016, 02:46 AM   #4
wdecoster
Member
 
Location: Antwerp, Belgium

Join Date: Oct 2015
Posts: 97
Default

Perhaps a quick and dirty assembly with miniasm can give you an idea? https://github.com/lh3/miniasm
wdecoster is offline   Reply With Quote
Old 07-15-2016, 06:51 AM   #5
jpummil
Member
 
Location: Fayetteville, AR

Join Date: Apr 2014
Posts: 85
Default

Quote:
Originally Posted by wdecoster View Post
Perhaps a quick and dirty assembly with miniasm can give you an idea? https://github.com/lh3/miniasm
Thanks for the suggestion wdecoster I think I've avoided miniasm thus far because it appears to only output .gfa files? Kind of limits further evaluation of the assembly as most common tools seem to still only take .fasta.

Update: Found a note in another thread about converting .gfa to .fasta Trying it now...

awk '/^S/{print ">"$2"\n"$3}' in.gfa | fold > out.fa

Last edited by jpummil; 07-15-2016 at 09:41 AM. Reason: Additional info
jpummil is offline   Reply With Quote
Old 10-03-2016, 05:27 AM   #6
jsoghigian
Junior Member
 
Location: New Haven, CT

Join Date: Sep 2016
Posts: 1
Default

Quote:
Originally Posted by jpummil View Post
Thanks for the suggestion wdecoster I think I've avoided miniasm thus far because it appears to only output .gfa files? Kind of limits further evaluation of the assembly as most common tools seem to still only take .fasta.

Update: Found a note in another thread about converting .gfa to .fasta Trying it now...

awk '/^S/{print ">"$2"\n"$3}' in.gfa | fold > out.fa
About to try this method myself - jpummil, were you successful in estimating genome size from your raw reads?
jsoghigian is offline   Reply With Quote
Old 10-06-2016, 09:44 AM   #7
jpummil
Member
 
Location: Fayetteville, AR

Join Date: Apr 2014
Posts: 85
Default

Quote:
Originally Posted by jsoghigian View Post
About to try this method myself - jpummil, were you successful in estimating genome size from your raw reads?
The assembly itself using miniasm and the conversion script from gfa to fasta worked fine, though the assembly isn't as "good" as from Canu.

Still no really good way to estimate a genome size from the PacBio reads. Schatz put together a really nice tool called GenomeScope, but currently only works with Illumina reads.
jpummil is offline   Reply With Quote
Old 10-10-2016, 05:00 AM   #8
Markiyan
Senior Member
 
Location: Cambridge

Join Date: Sep 2010
Posts: 115
Lightbulb Try analysing the reads from the short inserts (multipass ones).

You can try extracting the long raw reads from the short library inserts, which pass the insert multiple times (CCS-like reads), doing self error correction, and than using kmer counter software designed for Illumina, 454 or Sanger data.

Also please be aware, that you may have to screen out the high copy number DNA (mitochondrial/plastid genomes) before doing kmer counting.

Also you may get some PCR-Free miseq data to complement your pacbio assembly. (Can be cheaper if your coverage is still too low).
Markiyan is offline   Reply With Quote
Old 10-06-2019, 08:15 AM   #9
kartika
Junior Member
 
Location: london

Join Date: Oct 2019
Posts: 1
Default

thanks you
kartika 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 01:16 AM.


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