Peak in transcript coverage profile

[TomHarrop]

Hi all,

I'm analysing some RNASeq results using the tophat–htseq–DESeq2 pipeline. The libraries were prepared using polyA selection and some of our samples had low RINs, so I was checking the data for evidence of 3′ bias. I used Qualimap 1.0's RNA-Seq QC tool to generate transcript coverage profiles for each of the samples, and bias is indeed evident in the samples that had low RINs. However, ALL the samples, including the ones with acceptable RINs, have a noticeable double peak in the coverage profile at 62.5–67.5 %.

I'm attaching a coverage plot from a sample with RIN of 7.5, one of the good ones. I'm curious whether anyone can explain the double peak or has observed it before.

Thanks,

Tom

PS the double peak looks a bit like batman to me...

[TomHarrop]

A quick update. I'm still stumped by this issue. I did try re-mapping the reads with a different GTF file and genome fasta file, both of which I have used in previous analyses without problems. The strange double-peak is still present. This seems to rule out a problem with the annotation.

The libraries were generated from very small amounts of RNA (picogram range, I believe) after amplification with the Ovation RNA Amplification System v2. The libraries themselves were prepared with the Ovation RNA-Seq System v2. I'm wondering if the first amplification step is somehow responsible for this uneven coverage profile. Has anyone else had any experience with that particular kit?

Thanks again,

Tom

[TomHarrop]

Problem solved!

Not sure if anyone's reading this but I thought I'd do the right thing and post the answer in case someone googles it up in the future.

Right at the start of chromosome 9 of the MSU7 Oryza sativa annotation, there are two annotated "genes" (LOC_Os09g00999 and LOC_Os09g01000) that overlap the 17S and 25S rRNA subunit genes. This annotation caused millions of rRNA reads (which weren't completely removed by PolyA selection) to map to these two genes. Normally, the GTF would only contain mRNA sequences so this wouldn't occur.

There's a bit more info including a nice diagram (Fig. 3) in this publication: http://www.sciencedirect.com/science...88754313002024. For my analysis I will mask the two genes to prevent problems downstream e.g. incorrect estimation of size factors by DESeq2 etc.