Hi I am runing into very frustrating issue. I have RNA seq pair end samples that i am processing with Rsubread package. For the sample nr1 test everyhting worked fine and i obtained feature counts which look reasonable. However, after alignment of sample nr 2, attempt to perform feature counting returns this message

nput files : 1 BAM file || || || || OV2_Sample_aligned.bam || || || || Paired-end : yes || || Count read pairs : yes || || Annotation : gencode.v38.annotation.gtf (SAF) || || Dir for temp files : . || || Threads : 1 || || Level : meta-feature level || || Multimapping reads : counted || || Multi-overlapping reads : not counted || || Min overlapping bases : 1 || || || \============================================================================//

i dont have a slightest idea why the program says the file is SAF ( it is GTF from encode and it worked for the first sample)

i used relatively straight forward code

''library(Rsubread)

bam_file <- "path/aligned_reads.bam" gtf_file <- "path/to/annotation.gtf"

counts <- featureCounts(files = bam_file, annot.ext = gtf_file, isPairedEnd = TRUE)"

what could be wrong???? GTF file is absolutely fine i did check it in text editor.

Any suggestions?

The help page for featureCounts says that the external annotation is by default assumed to be in SAF format. If annot.ext is in GTF format, then you need to set isGTFAnnotationFile=TRUE.

The default settings for each argument can be seen from ?featureCounts or by typing:

> args(featureCounts)
function (files, annot.inbuilt = "mm39", annot.ext = NULL, isGTFAnnotationFile = FALSE, 
    GTF.featureType = "exon", GTF.attrType = "gene_id", GTF.attrType.extra = NULL, 
    chrAliases = NULL, useMetaFeatures = TRUE, allowMultiOverlap = FALSE, 
    minOverlap = 1, fracOverlap = 0, fracOverlapFeature = 0, 
    largestOverlap = FALSE, nonOverlap = NULL, nonOverlapFeature = NULL, 
    readShiftType = "upstream", readShiftSize = 0, readExtension5 = 0, 
    readExtension3 = 0, read2pos = NULL, countMultiMappingReads = TRUE, 
    fraction = FALSE, isLongRead = FALSE, minMQS = 0, splitOnly = FALSE, 
    nonSplitOnly = FALSE, primaryOnly = FALSE, ignoreDup = FALSE, 
    strandSpecific = 0, juncCounts = FALSE, genome = NULL, isPairedEnd = FALSE, 
    countReadPairs = TRUE, requireBothEndsMapped = FALSE, checkFragLength = FALSE, 
    minFragLength = 50, maxFragLength = 600, countChimericFragments = TRUE, 
    autosort = TRUE, nthreads = 1, byReadGroup = FALSE, reportReads = NULL, 
    reportReadsPath = NULL, maxMOp = 10, tmpDir = ".", verbose = FALSE)

Thank you so much! this helped :) I a am aware that all of it is written somewhere but sometimes very obvious things are hidden in plainsight. Thanks