I am currently analysing a RNAseq dataset of human lymphocytes containing a vector (which has mCHERRY as reporter gene and a different gene A).
In order to include the annotation for the new genes I have followed the instructions of Alexander Dobin (https://groups.google.com/forum/#!topic/rna-star/FGQRotrCB1Q) and it seems like STAR (v2.7.3a) is detecting them and including them in the BAM files.
However, when I try summarising the BAM files using featureCounts (I am using Rsubread 1.34.4 in my university's computer cluster and subread 2.0.0 at my local machine) these genes are not listed (not even with 0 counts). I have also tried adding an additional line for each gene changing "exon" for "gene", but it does not solve the issue.
Do you have any suggestions of why this might be?
This is the code I am using in the Rsubread:
feature_counts <- featureCounts( nthreads = 6, isGTFAnnotationFile = TRUE, GTF.featureType = "exon", annot.ext = "genome/gencode.v31.primary_assembly.annotation.gtf", GTF.attrType = "gene_id", primaryOnly = TRUE, isPairedEnd = TRUE, files = paste( "alignment/", list.files( path = "alignment/", pattern = ".bam$", recursive = TRUE, include.dirs = TRUE )[-1], sep = "" ) )
And in subread:
featureCounts -T 2 -t exon -a ./gencode.v31.primary_assembly.annotation.gtf -t exon -g gene_id -M --primary -p -o ./count_matrix/counts.txt \./alignment_test/AG*/*.bam