Dear Gilbert I'm confused you talk about the absence of an easy mode. However, in the web site of EvidentialGene there is a script based approach which looks as easy as:
evigene/scripts/prot/tr2aacds.pl -mrnaseq Trinity.fasta
It works for me and I obtained a transcriptome reduction of a 90% (from 1.2M transcripts to more or less 120k). I have the folders with the okay (and okalt) and dropset datasets. All looks fine.
#t2ac: EvidentialGene tr2aacds.pl VERSION 2013.07.27
# Class Table for ../../Trinity.trclass
class okay drop okay drop
althi 1.8 9.1 18541 91851
althi1 0.5 3.7 5585 37501
althia2 0 1 0 10354
altmfrag 0 0.2 803 2626
altmfraga2 0 0 112 191
altmid 0.4 1.7 4930 17350
altmida2 0 0.1 498 1296
main 4.2 7.1 42273 72170
maina2 0 0.2 752 2577
noclass 8.1 55.7 81985 560451
noclassa2 0 0.1 443 1521
parthi 0 4.2 0 42710
parthi1 0 0.6 0 6316
parthia2 0 0.1 0 1915
---------------------------------------------
total 15.5 84.4 155922 848829
=============================================
# AA-quality for okay set of ../../Trinity.aa.qual (no okalt): all and longest 1000 summary
okay.top n=1000; average=2037; median=1740; min,max=1363,9064; sum=2037028; gaps=0,0
okay.all n=125453; average=188; median=112; min,max=41,9064; sum=23667671; gaps=0,0
I have missed something? I should tweak something? I think my only weakness is to use only Trinity and not add other assembler as you suggested.
I'm using this "curated" transcriptome to my downstream analysis. If it is a wrong procedure plz correct me I'm learning and I appreciate all feedback to improve my analysis.
-- Pablo
For de novo assembly, I am currently in the process of running Velvet/Oases as an example of de novo assembly (takes a long time) and was considering running IDBA-Tran and SOAPdenovo-Trans also. For these I am only supplying my raw fastq reads. I am unable to use Trinity due to lack of linux or a server and because Galaxy Trinity doesn't seem to work properly.
I recommend you this paper: Informatics for RNA Sequencing: A Web Resource for Analysis on the Cloud http://journals.plos.org/ploscompbiol/article?id=10.1371/journal.pcbi.1004393
In relation to your questions:
I work with RNA-seq Data and have never used these approaches. working with linux in a server really makes everything easier. I have aways aligned with STAR or Trinity (De Novo), quantified using Kallisto/RSEM and analysed using Bioconductor tools (R).
I have see people using velvet for genome but never saw for transcriptomics. And I am not a fan of Cufflinks.
Hi Tiago,
Thank you for the reply. I took a look at the paper you referenced but it seems to be a little vague/skim over the bit that I'm most interested in. There also seems to be a lot of emphasis on Cufflinks in there. I'm actually at this stage not too concerned about the quantification of the different transcripts, at this stage I'd just like to visualize the data in a genome browser (as i have just a couple of genes that we're really interested in to focus in on). Just to clarify, I assume you use STAR for genome-guided and Trinity for De novo assembly of transcripts? And then do you combine the results of the 2 somehow to get a final list of transcripts that you have confidence are real? With your approach, do you get out of the end a set of transcripts that are known and also some that are novel? Presumably you have to take your De novo transcripts and map them back on to known transcripts somehow?
I would recommend IGV or Golden Helix Genome Browse, you can load your BAM file in both, the problm with golden helix is that you cant add your own reference, so if you do not work with mouse or human is better use IGV. Using the UCSC genome browse, the best way is to producing bigWig files, put it in a http server and load the tracks (this can be little overwhelming). There is also a good R package for visualization named gviz.
I only used STAR for mouse since the genome is very well annotated. I use trinity for some neglected disease without good genome reference.
If you want to map the trinity output back into your reference you can use BLAT.
hi,
I use StringTie without the
-eoption to enable known as well as novel transcript assembly. And it works well in picking out both. At least the known ones I can confirm. And for a cell line we had knocked out a gene by CRISPR, StringTie could assemble the affected transcript isoform. The cmd and the assembled isoforms -(input is STAR aligned)
Top blue is the canonical form, green is the control cell line and the next two are effect of knock-down
Hi Amitm,
Thank you for the reply. Do you have any experience using the newer version of Stringtie? It's just I seem to have used exactly the same parameters/commands as you and it gives me that output you see in the image above. I don't know why I'm not getting results like you have shown. I'm confused!
hi,
I have used the ver as noted in the code and older as well. All worked well. Haven't used newer ver. Maybe you could check the BAM file and see if you see spliced-alignment in IGV. A RNA-seq file aligned by a splice aware aligner, when loaded into IGV would show islands of reads linked by connecting bridges between the exons for any given gene. If you see this in your BAM on IGV then the transcript assembler should work. Also if you choose the Sashimi plot option in IGV you should see connecting loops (indicating split alignment).
Hi Amitm,
So you hit the nail on the head. There was a lack of communication in our lab. I thought that my reads had been aligned (by my colleague) with Subread using the "subjunc" function. Turns out it hadn't and therefore didn't have exon junction info. Once aligned properly and checked as you suggested via IGV, it ran perfectly with stringtie and I can see the results I was after. Feel like a bit of a dunce but thank you for helping me get to the bottom of my problem!
hi Amitm, kindly check this error.i did not get it and solve it .
./stringtie G1_sorted.bam -B -o G1.gtf -G Triticum_aestivum.IWGSC.42.gtf -p 4 -C G1.refs.gtf -A G1.abund.tab -WARNING: no reference transcripts were found for the genomic sequences where reads were mapped! Please make sure the -G annotation file uses the same naming convention for the genome sequences
Kindly set this command for my sample.
as it is written in the error message, check if -G annotation file uses the same naming convention for the genome sequences. If your bam files are aligned to a genome with
chrprefix then your GTF should also havechr, it usually happens when the reads are aligned to a genome fasta from UCSC and the GTF is from ENSEMBL.