Finding assembled transcript sequences from StringTie output
0
1
Entering edit mode
8 months ago
Jen ▴ 20

Hi All, I am currently using StringTie to identify novel transcripts (MSTRG ID transcripts) from RNAseq data in mouse. I have the gene.abund.tab output files from StringTie which gives the following info: Gene ID, Gene Name, Reference, Strand, Start, End, Coverage, FPKM, TPM. So from this file, I have the genomic loci where the assembled transcript maps to, but I don't know where to find the actual transcript sequence.

Also, an additional thing that is confusing me is that when I search the genomic location in the mouse genome using IGV, I'm finding that MSTRG transcripts are mapping to the full length known genes. How do I tell what is novel about them compared to the known annotated transcript? I'm guessing that if I know the actual sequence I can find this out.

-Jen

ADD COMMENT
0
Entering edit mode

Hi Jen! It would be helpful if you posted the commands you used for building the transcriptome. I'm assuming you assembled your sample in separate and then you merged them (using the reference annotation) using stringtie --merge. You are correct in assuming MSTRG are novel transcripts. However, that does not mean that they are novel loci. It is unclear to me if you are searching for novel genes in the transcriptome or you just want novel transcripts that might be in already known genes.

If you want to get the sequences you can run gffread on your GTF output file

ADD REPLY
0
Entering edit mode

Sorry, I should've been much more specific. I mapped the reads using STAR and the ran the following:

stringtie-2.1.4/stringtie -p 40 -G genes/gencode.vM25.primary_assembly.annotation.gtf -o Sample1.gtf star/Sample1_Aligned.sortedByCoord.out.bam


I did the above for all samples, then ran the following:

stringtie-2.1.4/stringtie --merge -p 40 -G genes/gencode.vM25.primary_assembly.annotation.gtf -o All_Samples.gtf mergelist.txt
stringtie-2.1.4/stringtie -e -B -p 40 -A Sample1_gene_abund.txt -G All_Samples.gtf -o ballgown/Sample1/Sample1.gtf Sample1.bam


I am currently looking at those MSTRG IDs that have no associated Gene ID, thus should be novel transcripts. From what I understand (I could be wrong), they could have no Gene ID because 1) they are new splice variants of a known gene (I would think these would have a Gene ID, so this may not be a valid reason) 2) They are a combination of 2 or more known known loci 3) They map to a completely new loci

ADD REPLY
0
Entering edit mode

I went and looked into one of my old GTF files and what I gathered was that your locus will only have the
"gencode.vM25.primary_assembly.annotation.gtf" gene ID if that locus does not contain new transcripts from your data.

One easy way to check for new loci would be running gffcompare:

gffcompare your_merged_transcriptome.gtf -o merged_transcriptome -r referece_gene_models.gtf


This will output a file with your_merged_transcriptome.gtf with an additional attribute in the last collum (gene_name), if there is an intersecting reference gene. Hope this helps

ADD REPLY
0
Entering edit mode

I'm a little confused about what files to use as input to run the gffread on my output .gtf file. I'm looking at the GFF utilities page at the code to "Extract transcripts sequences". Not sure if this is even the correct code to run....

gffread -w transcripts.fa -g /path/to/genome.fa transcripts.gtf

If it is, which files would I use? It seems that genome.fa would probably be the mouse genome file I used for the alignment, but unsure of the rest.

ADD REPLY
0
Entering edit mode

Your comand look fine to me. You have to use the same mouse genome that was used for the alignment and transcriptome assembly. What that command does is basically looking at your GTF file coordinates and extract the corresponding bases from the genome file

ADD REPLY
0
Entering edit mode

-w transcripts.fa <- is this the output file from gffread containing the transcript sequences I'm looking for? -g /path/to/genome.fa = I know this is the genome file used for alignment and assembly transcripts.gtf <- is this the output .gtf from StringTie

Sorry is this is obvious....

ADD REPLY
0
Entering edit mode

Yes you are correct in that assumption

ADD REPLY

Login before adding your answer.

Traffic: 2665 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6