SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
BFAST Web Server for the alignment of ABI SOLiD Data nilshomer Bioinformatics 2 05-18-2010 09:17 AM
ABI SOLiD data filtering and conversion to base-space PRJ SOLiD 5 12-15-2009 06:55 AM
ABI SOLiD $10k Grant Program: 60GB Mappable Sequencing Data ECO SOLiD 3 09-08-2009 05:53 PM
Abi Solid WT Sample Data vs. CLC Bio Genetics Workbench 3.0.0 Blaize RNA Sequencing 2 02-23-2009 04:44 AM
ZOOM released (supporting both Illumina data and ABI SOLiD data) spirit Bioinformatics 2 08-21-2008 07:48 AM

Reply
 
Thread Tools
Old 11-07-2010, 06:39 PM   #1
repinementer
Member
 
Location: asia

Join Date: Dec 2009
Posts: 80
Default ABI-SOLID data with Bowtie-0.12.7 and TopHat-1.1.2

Hi guys.
I successfully ran Solid data (50 bp reads) with current bowtie and tophat versions. But surprisingly it gave the following stats. especially only 21% aligning results!. Any one have any idea why aligning failed miserably ?

The command I used

Code:
tophat -o RHE014_Tophat_Output -C --library-type fr-secondstrand --segment-length 50 /home/choijk/Software/bowtie-0.12.7/indexes/hg18_c SL001_R00089_RHE014_01pgx2_F3_3_1
Following is the tophat log file description

Code:
$ cat filezqaMDd.log 
# reads processed: 97340093
# reads with at least one reported alignment: 20468407 (21.03%)
# reads that failed to align: 76790363 (78.89%)
# reads with alignments suppressed due to -m: 81323 (0.08%)
Reported 34579945 alignments to 1 output stream(s)
repinementer is offline   Reply With Quote
Old 11-07-2010, 09:29 PM   #2
KevinLam
Senior Member
 
Location: SEA

Join Date: Nov 2009
Posts: 203
Default

if you can maybe try using Bioscope/mapreads for the mapping.
it improves the odds.
~50% is expected for Whole transcriptome.
So you are getting lower than average..
KevinLam is offline   Reply With Quote
Old 11-07-2010, 09:52 PM   #3
repinementer
Member
 
Location: asia

Join Date: Dec 2009
Posts: 80
Default

Hmmm I already have results with Bioscope but I want to run these by using Tophat with out giving any known refseq annotation in order to find new transcripts.

Anyone suggestions on Tophat/bowtie
repinementer is offline   Reply With Quote
Old 11-07-2010, 09:57 PM   #4
KevinLam
Senior Member
 
Location: SEA

Join Date: Nov 2009
Posts: 203
Default

I see. How's the mapping % with bioscope then?
KevinLam is offline   Reply With Quote
Old 11-07-2010, 10:07 PM   #5
repinementer
Member
 
Location: asia

Join Date: Dec 2009
Posts: 80
Default

is around ~ 80 %
repinementer is offline   Reply With Quote
Old 11-08-2010, 12:12 PM   #6
adarob
Member
 
Location: Berkeley, CA

Join Date: Jul 2010
Posts: 71
Default

Can you try using fr-firststrand?
adarob is offline   Reply With Quote
Old 11-08-2010, 07:07 PM   #7
repinementer
Member
 
Location: asia

Join Date: Dec 2009
Posts: 80
Default

but fr-secondstarand is meant for SOLID data right ?
repinementer is offline   Reply With Quote
Old 11-09-2010, 12:08 PM   #8
adarob
Member
 
Location: Berkeley, CA

Join Date: Jul 2010
Posts: 71
Default

In general yes, but it depends on the protocol you used.
adarob is offline   Reply With Quote
Old 11-09-2010, 01:04 PM   #9
bacdirector
Member
 
Location: Pittsburgh, PA

Join Date: Oct 2009
Posts: 12
Default

if you're using 50b reads, try trimming them to 35b in colorspace (i.e., use only the first 35 bases). Most of the -1's (no-calls) and sequencing errors occur between bases 36 and 50.

We've DOUBLED our mapping using bioscope this way (from 25 million read to 50million reads).

There is a separate thread on this topic, but I'd be happy to hear back from you if you try this.

Let me know if you'd like to peek at our source code for trimming in colorspace.
bacdirector is offline   Reply With Quote
Old 11-09-2010, 05:11 PM   #10
rdeborja
Member
 
Location: Toronto

Join Date: Aug 2008
Posts: 42
Default

I'm curious to know if people are quoting the mapping percent as identified in the alignment.txt report or if they are calculating it based on total mapped reads as a percent of total reads. The Bioscope alignment report shows mapped reads as 100% then all figures others based on this.

I've been reporting total mapped reads, unique aligned reads, ribosomal and unmapped reads as a percent of total reads. Total mapped is usually 60-70%.

With the trimming to 35bp, is it necessary since Bioscope is doing a seed and extend on the reads. I find my average aligned read length to be ~40bp with a size frequency plot being bimodal at 25bp and 50bp.
rdeborja is offline   Reply With Quote
Old 11-09-2010, 06:01 PM   #11
repinementer
Member
 
Location: asia

Join Date: Dec 2009
Posts: 80
Default

@bacdirector: Yes I would be happy to do that. COuld you please provide me the code. My SOLID .csfasta data format looks like this
Quote:
# Wed Mar 10 00:25:29 2010 /share/apps/corona/bin/filter_fasta.pl --output=/data/results/sl001/SL001_R00089/RHE012_01pgx2/results.F1B1/primary.20100310065316729 --name=SL001_R00089_RHE012_01pgx2 --tag=F3 --minlength=50 --mincalls=25 --prefix=T /data/results/sl001/SL001_R00089/RHE012_01pgx2/jobs/postPrimerSetPrimary.2745/rawseq
# Cwd: /home/pipeline
# Title: SL001_R00089_RHE012_01pgx2
>853_10_97_F3
T.....023..2.10..120.3.2010.031...2.1.30.22001..00.
>853_10_111_F3
T.....113....33..003.0.010..100...2.0..2.03002..02.
>853_10_157_F3
T.....230..2.00..330.2.1313.231...3.1.10.02031..10.
>853_10_194_F3
T.....031....32..323.0.322..100...3.0..1.30313..23.
>853_10_221_F3
repinementer is offline   Reply With Quote
Old 11-24-2010, 02:30 AM   #12
plabaj
Member
 
Location: Vienna

Join Date: Oct 2010
Posts: 86
Default Problem with TopHat and ABI SOLiD

I don't know if you are aware but the current version TopHat is using different algorithm than was described in TopHat paper from 2009. The current algorithm is described in supplement to Cufflinks paper (Trapnell 2010).

The most important change is that read is split to segments, and discovering of the splice junctions is based where these segments aligned. New TopHat is optimsed for >=75bp reads, in this case each read is divided to 3 segments each 25bp.

You run TopHat with setting the segment length to 50 (--segment-length 50), which means that there will be just ONE segment, thus such setting cannot discover any splice junction.

And here is the problem with current TopHat, it seems to be NOT designed for ABI SOLiD reads. You have two options for 50bp reads:
- use default settings, and be aware that not all splice junctions will be discovered
- set --segment-length 16, but then 16bp segments will align everywhere and in many cases will be discarded, so again many splice junctions will not be discovered and many false positives will be found

For 36bp reads situations is even worse.
Old versions of TopHat supported such short reads but didn't support color space reads

In my opinion, so far, there is no proper software for discovering new transcripts, or even assembling properly existing ones, for ABI SOLiD data. If you know any please let me know.
plabaj is offline   Reply With Quote
Old 12-15-2010, 08:30 PM   #13
tsucheta
Member
 
Location: Arlington, TX

Join Date: Nov 2009
Posts: 17
Default

@plabaj:
Very True!! My experience with tophat has been disaster so far.. While bowtie gives way more alignments; using tophat for single end reads, it reports only 24% alignment. My reads size is 50 bases. The junctions.bed file produced is extremely unreliable since it reports only 851 junctions in < 10% of the scaffolds. The remaining is unreported! On the top of it, most of the junctions are overlapping. Is there a possible reason where it gets things wrong?
tsucheta is offline   Reply With Quote
Old 12-16-2010, 12:02 AM   #14
plabaj
Member
 
Location: Vienna

Join Date: Oct 2010
Posts: 86
Default

@tsucheta

As I wrote before.
New algorithm for finding splice junctions implemented in TopHat is responsible for that.
It starts working properly for reads of length 75bp.

If your reads are not ABI SOLiD try to install older version of TopHat (you have to find out which one analysing the versions changes).
Maybe it will help.

Or if you are interested just in transcripts expresion use BowTie against transcripts sequences (not genomic sequence like for TopHat)
__________________
Pawel Labaj
plabaj is offline   Reply With Quote
Old 12-16-2010, 01:19 AM   #15
xinwu
Member
 
Location: Beijing

Join Date: Jul 2010
Posts: 33
Default

Quote:
Originally Posted by tsucheta View Post
@plabaj:
Very True!! My experience with tophat has been disaster so far.. While bowtie gives way more alignments; using tophat for single end reads, it reports only 24% alignment. My reads size is 50 bases. The junctions.bed file produced is extremely unreliable since it reports only 851 junctions in < 10% of the scaffolds. The remaining is unreported! On the top of it, most of the junctions are overlapping. Is there a possible reason where it gets things wrong?
What's your segment size?
xinwu is offline   Reply With Quote
Old 01-29-2011, 09:44 PM   #16
tsucheta
Member
 
Location: Arlington, TX

Join Date: Nov 2009
Posts: 17
Default

@Xinwu, My segment length was default 25.

Last edited by tsucheta; 01-30-2011 at 01:31 PM.
tsucheta is offline   Reply With Quote
Old 03-31-2011, 11:08 AM   #17
tonge
Junior Member
 
Location: canada

Join Date: Mar 2011
Posts: 5
Default

@plabaj (and anyone else),
I have read your comment with interest and some concern. I plan to generate SOLID RNA-seq data (paired end 50bp+25bp) and wonder whether Tophat is any better at handling paired end SOLID data. I was planning to utilise Tophat for helping with new transcript/exon discovery. BTW this is all with mouse, so have annotated genome to work with.
Thanks.
tonge is offline   Reply With Quote
Old 04-19-2011, 09:10 AM   #18
mackan
Junior Member
 
Location: Sweden

Join Date: Dec 2009
Posts: 4
Default

Hi All,

I am trying to align paired-end Solid whole-transcriptome reads (50+35 FR) using Bowtie 0.12.7.

Strangely, when mapping paired-end, no read pairs (other than a few repeat-regions) map to the genome. In addition, no read pairs map to the transcriptome. (Ensembl genes-based reference).

When mapping individual reads, about 30% of 50bp reads and 20% of 35bp reads map successfully.

I must be doing something wrong. Using --ff changes nothing (reads are actually FR. The older Solid mate-pairs were FF) Read csfasta files are 100% pair matched (with read name suffixes _F3 and _F5-BC).

Examples (using only 1k reads but the results hold for the full set - as well as with exon-only Ensemble genes as reference, one fasta element per gene)
PE:
bowtie -f -p 8 --fr -C -S --sam-nohead --sam-nosq -s 1000000 -u 1000 --Q1
$dir/F3_QV.qual --Q2 $dir/F5_QV.qual -1 $dir/F3.csfasta -2
$dir/F5.csfasta ~/p2/indexes/bowtie/hg19_c align.sam
# reads processed: 1000
# reads with at least one reported alignment: 1 (0.10%)
# reads that failed to align: 999 (99.90%)
Reported 1 paired-end alignments to 1 output stream(s)

The aligned pair:
17_213_1598_F3 67 chr2 154876116 255 48M = 154876139 56 CACACACACACACACACACACACACACACACACACACACACACACACA LbRL\TMUVBH\TQ.8[[bOM```LG]^_LJ^\_SN]bcca_aaabbW XA:i:1 MD:Z:48 NM:i:0 CM:i:1
17_213_1598_F5-BC 131 chr2 154876140 255 33M = 154876115 -58 CACACACACACACACACACACACACACAACAGC @NcPIAE`c^E:S_^KI]`cPGZZSED)!E1!3 XA:i:0 MD:Z:29A3 NM:i:1 CM:i:5

F3-ends (50b) only
[markus@q34 run]$ bowtie -f -p 8 -C -S --sam-nohead --sam-nosq -s
1000000 -u 1000 -Q $dir/F3_QV.qual ~/p2/indexes/bowtie/hg19_c
$dir/F3.csfasta algn.sam
# reads processed: 1000
# reads with at least one reported alignment: 319 (31.90%)
# reads that failed to align: 681 (68.10%)
Reported 319 alignments to 1 output stream(s)

F5-ends (35b) only
[markus@q34 run]$ bowtie -f -p 8 -C -S --sam-nohead --sam-nosq -s
1000000 -u 1000 -Q $dir/F5_QV.qual ~/p2/indexes/bowtie/hg19_c
$dir/F5.csfasta algn.sam
# reads processed: 1000
# reads with at least one reported alignment: 287 (28.70%)
# reads that failed to align: 713 (71.30%)
Reported 287 alignments to 1 output stream(s)

Running tophat generated a decent result, so plenty of reads should map in pairs.

I'd be grateful if anyone can help me out, or provide any hint.

Edit: I got it working now, I might have forgotten about the --fr for the full data set.

Last edited by mackan; 04-20-2011 at 04:46 AM.
mackan 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:12 PM.


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