RNAseq reads mapped to the gene
2
0
Entering edit mode
@myprogramming2016-9741
Last seen 6.8 years ago

Hi,

I work on plant species. I have used STAR to map RNAseq reads and featureCounts to get expression values. I would like to counts the number of reads map to the genes and outside of genes. Is there any tool or script to get these estimates?

Thanks

star featurecounts edgeR RNAseq • 2.1k views
ADD COMMENT
0
Entering edit mode

Did I understood correctly that you want 2 numbers: the number of reads mapping to a gene feature and those which don't?

ADD REPLY
0
Entering edit mode

Yeah, that's true.  I would like to get this estimate for each gene and each of my Individual.

Thanks

ADD REPLY
1
Entering edit mode

So you have the total number of reads which you mapped, and you have the counts generated by featureCounts. I guess you can just make the sum of those assigned to a gene and consider those as counted, and subtract that from the total and take those as not-in-a-gene. But that's an oversimplification, since multimapping reads can be in a gene and will not be counted. Depends on your application what you want to do with those.

ADD REPLY
1
Entering edit mode

So now you want it per gene and per individual? So that's just what featureCounts does?!

ADD REPLY
0
Entering edit mode

Thanks!

 

ADD REPLY
2
Entering edit mode
phil.chapman ▴ 150
@philchapman-8324
Last seen 7.5 years ago
United Kingdom

Since featureCounts gives you a matrix of reads counts with genomic feature as the row name and sample as the column name, you can calculate the number of reads that have been mapped to the features defined in your annotation by calculating the sum of your matrix.

fc <- featureCounts(files=bam_files,
                    annot.ext='Homo_sapiens.GRCh37.75.gtf',
                    isGTFAnnotationFile=TRUE,
                    GTF.featureType="exon",
                    GTF.attrType="gene_id",
                    useMetaFeatures=TRUE,
                    allowMultiOverlap=FALSE,
                    nthreads=16,
                    strandSpecific=2,
                    countMultiMappingReads=FALSE)
reads_in_gtf_features <- sum(fc$counts)

It's well worth throughly reading the documentation for featureCounts since there are a lot of parameters and how you specify these will make a big difference to your end result. 

 

ADD COMMENT
0
Entering edit mode
@myprogramming2016-9741
Last seen 6.8 years ago

Thanks for your help. I will go through the featureCounts manual. 

I have used linux version of featureCounts.  

 

 

ADD COMMENT

Login before adding your answer.

Traffic: 841 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