featureCounts reports zero counts when using GTF Repeatmasker track from UCSC and Bowtie2 SAM mapped RNA-seq reads

Question

[Gema Sanz]

Hi,

I want to calculate the read counts mapping to repeats in my RNA-seq data. I mapped the reads using Bowtie2 to indices of repeats from UCSC Repeatmasker track in fasta format. Next, I have to calculate the read counts per repeat type, so I tried with featureCounts tool (in Galaxy server). But the output gives me zero counts for all the repeat types.

¿Any ideas? I thought maybe some IDs or columns are not matching but I don't know.

Thank you very much in advance, Gema

[h.mon]

Probably featureCounts is discarding multi-mapped reads, hence zero counts for repeats. You have to keep only one copy of the repeats on the genome, or count with RSEM or Salmon, which will try to estimate by expectation-maximization how to distribute the counts of multi-mapped reads.

[Gema Sanz]

Thanks for your answer h.mon! In Bowtie2 I used default mode (search for multiple alignments, report the best one), I'm not sure what you mean with "keep only one copy of the repeat". Reads are mapped directly to repeat indices not to human genome. I think there is an option in featureCounts to keep the multi-mapped reads, I will check.

Answer 1

[Gema Sanz]

I tried the -M option in featureCounts (allow multi-mapped) but still I'm getting zeros everywhere:

Geneid example.sam AluSp 0 AluY 0 L2b 0 L1PA10 0 L1PA2 0 L1MB7 0 ERVL-E-int 0

GTF repeatmasker file looks like this:

Seqname Source Feature Start End Score Strand Frame Attributes
chr1 hg19_rmsk exon 16777161 16777470 2147.000000 + . gene_id "AluSp"; transcript_id "AluSp";

And the sam file:

QNAME FLAG RNAME POS MAPQ CIGAR MRNM MPOS ISIZE SEQ QUAL OPT @HD VN:1.0 SO:unsorted
@SQ SN:hg19_rmsk_AluSp LN:310
@SQ SN:hg19_rmsk_AluY LN:289
@SQ SN:hg19_rmsk_L2b LN:1040
@SQ SN:hg19_rmsk_L1PA10 LN:2090

My feeling is that something is not matching but I don't know what...