Hi, I'm using STAR to align RNA-seq data to mm39. I am using --quantMode geneCounts as an option and the results I get are in the (very impractical for my purpose) 'ENSMUSG' format. Is there a way to get gene names instead of the Ensemble gene IDs?
Additionally, I'm confused as to why I get three columns of counts for each gene (as shown below) -- I'm only aligning one forward and one reverse fastq file, so shouldn't I be getting one set of reads per gene?
N_unmapped 16802670 16802670 16802670
N_multimapping 4055291 4055291 4055291
N_noFeature 22948357 58749681 25405703
N_ambiguous 3018494 37339 1274777
ENSMUSG00000102628 0 0 0
ENSMUSG00000100595 0 0 0
ENSMUSG00000097426 0 0 0
ENSMUSG00000104478 0 0 0
ENSMUSG00000104385 0 0 0
ENSMUSG00000086053 21 25 0
If it helps, the code I used to index the genome and align my fastq files is the following:
STAR --runMode genomeGenerate --genomeDir mm39index --genomeFastaFiles /path/to/file/Mus_musculus.GRCm39.dna.primary_assembly.fa --sjdbGTFfile /path/to/file/Mus_musculus.GRCm39.104.gtf --runThreadN 16
STAR --runThreadN 16 --genomeDir path/to/mm39index --readFilesIn blah_1.fastq.gz blah_2.fastq.gz --readFilesCommand zcat --outSAMtype BAM SortedByCoordinate --quantMode GeneCounts --outFileNamePrefix alignments/blah-alignment
Thanks in advance for the help!
2 answers
The three columns correspond to (I forget the exact order) stranded, reverse stranded, and unstranded library types.
It's recommended to keep them as gene ids and then merge a list of gene names because gene ids are unique, but gene names aren't. You can get a id to name mapping from the gtf file you used, or a database like biomart.
If all you care about is quantification Salmon is usually recommended these days. It gives more accurate quantification due to being able to quantify abundances at the transcript level, allowing corrections for things like transcript length bias.
HTseq could be another option for feature count rather than usisng quant mode in STAR. I tried that berfore but ended up using HTseq. The HTseq outputs can be combined to count matrix easily and feed to analysis pipeline.
Log in to answer this question.
-
STAR solo segmentation fault after 'started Solo counting'
written by bp22 •
-
No Unmapped Reads in bam file?
written by erik.burchard •
-
Weird output of geneCount for STAR?
written by FantasticAI •
-
Low mapping frequency on STAR
written by a_bis •
-
How should I deal with my RNAseq result with high % of mutiple aligned reads?
written by Simon Ahn •
-
How to deal with high % of reads mapped to multiple loci on STAR?
written by Simon Ahn •
-
Large Count Discrepancy of Key Gene Between STAR and HISAT2
written by dk0319 •
-
ReadsPerGene.out.tab column one contains gene ID AND gene name
written by lalina1606 •
-
STAR --quantMode GeneCounts always gives 0 counts
written by mgmohsen •
-
Counting No. of Reads per gene from RNA-Seq with STAR
written by skevalkumar •
