SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
Star with cufflinks Nino Bioinformatics 20 12-12-2016 09:06 AM
STAR + Cufflinks. Cufflinks hanging. BAM XS error? gogodidi Bioinformatics 3 12-20-2015 09:19 PM
using STAR+Cufflinks for transcript assembly turns unstranded RNA-seq to stranded? runxuan Bioinformatics 5 07-22-2015 03:24 AM
Reassign mapping quality in STAR / GATK/ Cufflinks compatibility cbaudo Bioinformatics 0 04-09-2015 06:50 AM
TimeLogicís latest FPGA technology to be coupled with Geneious Server Geneious Vendor Forum 0 01-24-2011 01:09 PM

Reply
 
Thread Tools
Old 02-08-2016, 01:06 PM   #1
harrike
Member
 
Location: St. Louis, MO

Join Date: Jun 2010
Posts: 29
Default STAR coupled with Cufflinks

I am working on a plant, in which the genome has been sequenced but the gene-annotation quality are not as good as I expected. I am trying to use STAR and Cufflinks to annotate novel transcripts from pair-end RNA-seq data. When I used the default settings, I got many alternative-spliced transcripts, some of which don't look like real. Therefore, for each novel splicing junction, I would like to have a filter, that is, at least ten reads (unique mapped or multi-mapped) support the junction. I don't know how to do that. What I am currently using is to do 2-pass mapping using STAR. After the 1st pass mapping, I will have the unsupported (with less than 10 reads) junctions removed from the *SJ.out.tab, and use this modified *SJ.out.tab as "annotated" junctions (--sjdbFileChrStartEnd value) for the 2nd pass mapping. Will this step improve the the mapping and get a more reliable list of novel transcripts.

Thanks!
harrike is offline   Reply With Quote
Old 02-09-2016, 03:26 PM   #2
alexdobin
Senior Member
 
Location: NY

Join Date: Feb 2009
Posts: 161
Default

Hi @harrike,

--outSJfilterCountUniqueMin 10 10 10 10 --outFilterType BySJout
will remove all junctions that are supported by <10 uniquely mapping reads. The 4 numbers are for different intron motifs: non-canonical, GT/AG, GC/AG, AT/AC.
You can also use --outSJfilterCountTotalMin 10 10 10 10 --outFilterType BySJout, which will remove junctions supported by <10 reads including unique and multimappers.

This option will work with 1-pass or 2-pass mapping.

There are also other splice junction filtering options - please check the manual for outSJfilter* otpions, they can all be used together in "AND" fashion.

Cheers
Alex
alexdobin is offline   Reply With Quote
Old 02-09-2016, 06:13 PM   #3
harrike
Member
 
Location: St. Louis, MO

Join Date: Jun 2010
Posts: 29
Default

Hi Alex,

Great ! Thank you very much for your reply. I am testing those settings.

I am curious about what happens to those uniquely mapped reads spanning a junction with <10 reads support. Are they considered as unmapped reads?

Thanks,

Rui
harrike is offline   Reply With Quote
Old 02-09-2016, 06:58 PM   #4
harrike
Member
 
Location: St. Louis, MO

Join Date: Jun 2010
Posts: 29
Default

Hi Alex,

I tried one of my libraries, using the "--outSJfilterCountTotalMin 10 10 10 10 --outFilterType BySJout". In the *SJ.out.tab file, some of the junctions are still of <10 reads support (see below). Should the sum of the column 7 and 8 be >=10?

scaffold2131 22287 22614 1 1 0 62 0 39
scaffold2131 22322 22614 1 1 0 3 0 23
scaffold2131 22504 22614 1 1 0 1 0 35
scaffold2131 22527 22614 1 1 0 1 0 16
scaffold2131 22702 22827 1 1 0 40 0 44
scaffold2131 22933 23197 1 1 0 31 0 44
scaffold2131 23273 23348 1 1 0 157 0 45
scaffold2131 23479 23582 1 1 0 94 61 44
scaffold2131 23675 23875 1 1 0 134 24 44
scaffold2131 23927 24382 1 1 0 150 0 44
scaffold2131 24460 24565 1 1 0 131 23 44
scaffold2131 24621 24694 1 1 0 2 0 34
scaffold2131 24621 24699 1 1 0 1 0 40

The command I used is:

STAR --genomeDir /home/harrikex/workspot/litchi/RNAseq/star_index/ --runThreadN 24 --readFilesIn 121224_I663_FCC1L8RACXX_L6_LITgqjTABRAAPEI-62_1.fq 121224_I663_FCC1L8RACXX_L6_LITgqjTABRAAPEI-62_2.fq --outSAMtype BAM Unsorted --outFilterMultimapNmax 20 --alignIntronMax 10000 --alignMatesGapMax 10000 --outSAMstrandField intronMotif --alignEndsType EndToEnd --outFilterIntronMotifs RemoveNoncanonical --outSJfilterCountTotalMin 10 10 10 10 --outFilterType BySJout --outFileNamePrefix FMO_121224_

Thanks,

Rui
harrike is offline   Reply With Quote
Old 02-10-2016, 01:02 AM   #5
Macspider
Member
 
Location: Vienna

Join Date: Feb 2016
Posts: 35
Default

The 2-pass pipeline for STAR works very fine, I used it many times. Moreover, STAR with Cufflinks scores among the best results in alignment-based quantification pipelines when dealing with not very specific samples, e.g. tissue samples or whole-organism samples.
Macspider is offline   Reply With Quote
Old 02-10-2016, 07:02 AM   #6
alexdobin
Senior Member
 
Location: NY

Join Date: Feb 2009
Posts: 161
Default

Quote:
Originally Posted by harrike View Post
Hi Alex,

I tried one of my libraries, using the "--outSJfilterCountTotalMin 10 10 10 10 --outFilterType BySJout". In the *SJ.out.tab file, some of the junctions are still of <10 reads support (see below). Should the sum of the column 7 and 8 be >=10?

...

Thanks,

Rui
Hi Rui,

sorry, I forgot that --outSJfilterCountUniqueMin and --outSJfilterCountTotalMin filters are used in the OR fashion, i.e. if one or the other is satisifed, the junction will be kept - so you have to specify both options. If you want only the Total filter, you need to make --outSJfilterCountUniqueMin numbers larger or -1, e.g. --outSJfilterCountUniqueMin -1 -1 -1 -1 --outSJfilterCountTotalMin 10 10 10 10 .

The reads that cross the junctions that do not pass this filter will be re-mapped forcing them not to cross these junctions. Most often this will result in soft-clipped reads, i.e. donor or acceptor part of the read over the prohibited junction will be dropped.

Cheers
Alex
alexdobin is offline   Reply With Quote
Old 02-10-2016, 07:14 PM   #7
harrike
Member
 
Location: St. Louis, MO

Join Date: Jun 2010
Posts: 29
Default

Quote:
Originally Posted by alexdobin View Post
Hi Rui,

sorry, I forgot that --outSJfilterCountUniqueMin and --outSJfilterCountTotalMin filters are used in the OR fashion, i.e. if one or the other is satisifed, the junction will be kept - so you have to specify both options. If you want only the Total filter, you need to make --outSJfilterCountUniqueMin numbers larger or -1, e.g. --outSJfilterCountUniqueMin -1 -1 -1 -1 --outSJfilterCountTotalMin 10 10 10 10 .

The reads that cross the junctions that do not pass this filter will be re-mapped forcing them not to cross these junctions. Most often this will result in soft-clipped reads, i.e. donor or acceptor part of the read over the prohibited junction will be dropped.

Cheers
Alex

Hi Alex,

Thanks for your explanation. I tried the combination of "--outSJfilterCountUniqueMin -1 -1 -1 -1 --outSJfilterCountTotalMin 10 10 10 10", it works. It has only half of the junctions (using the default setting) reported. It seems to have greater reliability.

One more question: After the 1st-pass mapping, using the above setting, only the junctions with >=10 reads support will be reported into the *SJ.out.tab file. As I have a couple of samples, I will get a few *SJ.out.tab files. In the 2nd-pass, all these *SJ.out.tab will be used as references (annotated junctions) for the mapping, do I still need to use the setting of "--outSJfilterCountUniqueMin -1 -1 -1 -1 --outSJfilterCountTotalMin 10 10 10 10" ? as it is possible, a few reads, <10 in one library, are mapped to a junction which is of >10 read support in another library; in this case, I still want those read not to be filtered out. On the other hand, if I don't use those settings, I feel many false positives will be introduced in the 2nd-pass mapping.

In other words, in the 2nd-pass mapping, will all the reads be preferentially mapped according the annotated junctions (provide in the *SJ.out.tab files using --sjdbFileChrStartEnd)? If this is true, I feel it is likely good without the >10 filter (--outSJfilterCountUniqueMin -1 -1 -1 -1 --outSJfilterCountTotalMin 10 10 10 10) in the 2nd-pass mapping.

Thanks,

Rui
harrike is offline   Reply With Quote
Old 02-11-2016, 07:57 AM   #8
alexdobin
Senior Member
 
Location: NY

Join Date: Feb 2009
Posts: 161
Default

Hi Rui,

the --outSJfilterCountUniqueMin -1 -1 -1 -1 --outSJfilterCountTotalMin 10 10 10 10 filter will not apply to the "annotated" junctions, i.e. junctions that are introduced via --sjdbFileChrStartEnd or --sjdbGTFfile options. So in the 2nd pass the junctions from your combined-SJ file will not be filtered out.
However, STAR will still be detecting the junctions which never made it to combined-SJ file, and you will have to use this filter to get rid of them.

Cheers
Alex
alexdobin is offline   Reply With Quote
Old 02-11-2016, 10:23 AM   #9
harrike
Member
 
Location: St. Louis, MO

Join Date: Jun 2010
Posts: 29
Default

Hi Alex,

Great, thanks for your clarification.

I decide to using the filter "--outSJfilterCountUniqueMin -1 -1 -1 -1 --outSJfilterCountTotalMin 10 10 10 10" for my 2nd-pass mapping.

Rui
harrike is offline   Reply With Quote
Old 02-25-2016, 04:38 PM   #10
biokloot
Junior Member
 
Location: SouthPacific

Join Date: Feb 2016
Posts: 2
Default

Dear Rui and Alex,
Many thanks for asking and explaining the SJout filter combinations so well!
I want to completely suppress detection of any and all SJ in the second pass. Do need to set --outSJfilterCountUniqueMin -1 -1 -1 -1 AND --outSJfilterCountTotalMin -1 -1 -1 -1 for the second passes to achieve that?
biokloot is offline   Reply With Quote
Old 02-25-2016, 05:25 PM   #11
harrike
Member
 
Location: St. Louis, MO

Join Date: Jun 2010
Posts: 29
Default

Hi biokloot,

To me, all "-1" does suppress all SJ, as "-1 means no output for that motif". But I am not sure how a read spanning an intron will be mapped in this case.

I am sure Alex will give a good clarification.

Rui
harrike is offline   Reply With Quote
Old 02-26-2016, 02:26 PM   #12
alexdobin
Senior Member
 
Location: NY

Join Date: Feb 2009
Posts: 161
Default

Quote:
Originally Posted by biokloot View Post
Dear Rui and Alex,
Many thanks for asking and explaining the SJout filter combinations so well!
I want to completely suppress detection of any and all SJ in the second pass. Do need to set --outSJfilterCountUniqueMin -1 -1 -1 -1 AND --outSJfilterCountTotalMin -1 -1 -1 -1 for the second passes to achieve that?
Hi @biokloot,

the best way to suppress detection of any new junctions in 2nd pass (i.e. junctions not detected in the 1st pass) is to use --alignSJoverhang 100000 (or any number > read length).
Alternatively, you can use
--outSJfilterCountUniqueMin -1 -1 -1 -1 --outSJfilterCountTotalMin -1 -1 -1 -1 --outFilterType BySJout
or
--outSJfilterOverhangMin -1 -1 -1 -1 --outFilterType BySJout
but it's slightly slower.

Cheers
Alex
alexdobin is offline   Reply With Quote
Old 02-28-2016, 12:08 AM   #13
biokloot
Junior Member
 
Location: SouthPacific

Join Date: Feb 2016
Posts: 2
Default

Thanks for the quick replies, Rui and Alex!
biokloot 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 09:27 AM.


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