SEQanswers

Go Back   SEQanswers > General



Similar Threads
Thread Thread Starter Forum Replies Last Post
How to align contigs without reference genome? yangfangisok Bioinformatics 5 06-17-2014 09:34 AM
how to align the contigs to the reference genome jjjscuedu Bioinformatics 1 06-05-2012 08:39 AM
Align primer against NGS contigs? -yl- Bioinformatics 3 11-20-2011 10:01 PM
Align reads to contigs ojy Bioinformatics 3 07-25-2011 09:16 AM
How to align contigs? azmicro 454 Pyrosequencing 10 10-06-2010 06:14 AM

Reply
 
Thread Tools
Old 02-08-2017, 01:18 PM   #1
sfh838t
Member
 
Location: Mountain Grove, MO, USA

Join Date: Apr 2014
Posts: 23
Default reads align but contigs do not??

Here is something that puzzles me:
I have 4000 reads that align to my reference genome. I can visualize them, they cover maybe 80 % or a little more of the reference sequence.
I assemble these reads into contigs, get something like 160 contigs. ok.
however, when I try to align these contigs, only 27 of them still match my reference.
huh? how is that possible?
Any thoughts?
sfh838t is offline   Reply With Quote
Old 02-08-2017, 01:55 PM   #2
GenoMax
Senior Member
 
Location: East Coast USA

Join Date: Feb 2008
Posts: 6,484
Default

Have you checked to see how much duplication you have in your data? 4000 data tends to have more duplicates if you don't have the right library insert size. Use clumpify.sh from BBMap to get an estimate of total duplicates and optical duplicates. You may not be getting a good assembly if there are too many duplicates.
GenoMax is offline   Reply With Quote
Old 02-08-2017, 05:59 PM   #3
Brian Bushnell
Super Moderator
 
Location: Walnut Creek, CA

Join Date: Jan 2014
Posts: 2,695
Default

Additionally, lengths and coverage are important here. Possibly, 99% of the bases map to the 27 contigs that match the reference, and the other contigs are short, degenerate contigs from errors or chimerism. Or maybe you're assembling various adapter sequences (another form of chimerism, in a way). It depends on how you're doing the alignment and assembly, and what kind of data you have. Can you describe the situation in more detail? Also, have you BLASTed the unaligned sequences, or compared them to known adapter sequences?
Brian Bushnell is offline   Reply With Quote
Old 02-09-2017, 05:10 AM   #4
sfh838t
Member
 
Location: Mountain Grove, MO, USA

Join Date: Apr 2014
Posts: 23
Default

these are single end reads from an RNA sequencing project. adapter sequences have been removed.
out of some 37 mil reads I get 4000 reads reads aligning to a just under 8000bp long virus that we know is present in the sample. But these are not siRNA in 20nt length range, these reads range from about 30 to 100 nt long.
I am sure there is duplication, though I do not know how to find or eliminate that and would appreciate any hints as to where to read about this or find tools to work with this.
out of about 160 contigs, 40 are very short, the same length as the shortest reads in fact.
the reads matching the virus (bwa aligned, IGV visualized) seem to visually cover most of the virus sequence, the 27 contigs cover maybe 30 % of the virus sequence.
Thanks for the help !
sfh838t is offline   Reply With Quote
Old 02-09-2017, 06:17 AM   #5
GenoMax
Senior Member
 
Location: East Coast USA

Join Date: Feb 2008
Posts: 6,484
Default

Quote:
Originally Posted by sfh838t View Post
I am sure there is duplication, though I do not know how to find or eliminate that and would appreciate any hints as to where to read about this or find tools to work with this.
@Brian's post gives directions on how to use clumpify.sh from BBMap to identify, mark/eliminate duplicates. They can be optical or all (PCR+ optical).
GenoMax is offline   Reply With Quote
Old 02-09-2017, 06:35 AM   #6
sfh838t
Member
 
Location: Mountain Grove, MO, USA

Join Date: Apr 2014
Posts: 23
Default

I did see Brians post, but did not check on that so far. Looks like I will have to take the plunge and install yet another toolbox. Let's see what will go wrong this time.....
sfh838t is offline   Reply With Quote
Old 02-09-2017, 06:40 AM   #7
GenoMax
Senior Member
 
Location: East Coast USA

Join Date: Feb 2008
Posts: 6,484
Default

Nothing should go wrong. As long as you have Java 1.7 or greater available nothing else is needed.
GenoMax is offline   Reply With Quote
Old 02-09-2017, 08:19 AM   #8
sfh838t
Member
 
Location: Mountain Grove, MO, USA

Join Date: Apr 2014
Posts: 23
Default

Ok, I did install and use this clumpify tool.
and it did not make any difference. it found some 300 duplicates, I assemble slightly fewer reads, to similar n50 and total nt count and of the contigs produced I had 1 less contig match than before.....
I have heard the term redundant vs non-redundant. just to make sure I understand: removing these duplicates would make my set of reads non-redundant? Or am I completely lost?
sfh838t is offline   Reply With Quote
Old 02-09-2017, 10:12 AM   #9
GenoMax
Senior Member
 
Location: East Coast USA

Join Date: Feb 2008
Posts: 6,484
Default

What setting did you use for dupedist? For a 4000 flowcell data you need to include
Code:
dedupe=t optical=t dupedist=2500 dupesubs=0 spantiles=f
to your clumpify command line. (dupesubs=0 will only look for perfectly matching reads without any errors).

If you did not use even one of those options then I suggest that you try again.
GenoMax is offline   Reply With Quote
Reply

Tags
assembly accuracy, contig, read

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 04:11 PM.


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