How does picard determine reads that map to correct strand?

[batman]

Hi,

I am curios to know how picards CollectRnaSeqMetrics tool determines the correct strand a read maps to. My guess is that it will look at the annotation supplied and determines the correctness of the strand based on that. Can anyone confirm?

[ak352]

According to http://picard.sourceforge.net/picard...#RnaSeqMetrics, the input is a SAM (Sequence Alignment Format) file.

This file contains the annotation related to on which strand and on which location on the reference a read is mapped to.
Refer to page 4 (0x10 bit of the flag field) of SAM format specification. I suppose the easiest way would be to just read this flag.

[batman]

Hi ak352,

Thanks for the reply. But I think I should have rephrased the question. What I needed to know was how can CollectRnaSeqMatrics know the correct strand when it calculates the field INCORRECT_STRAND_READS. i.e. how does it know the specific strand the original read came from, say when doing unstranded rnaseq mapping?

[kmcarr]

Hi ak352,

Thanks for the reply. But I think I should have rephrased the question. What I needed to know was how can CollectRnaSeqMatrics know the correct strand when it calculates the field INCORRECT_STRAND_READS. i.e. how does it know the specific strand the original read came from, say when doing unstranded rnaseq mapping?

The "(IN)CORRECT_STRAND_READS" is only meaningful and appropriate if your RNA-Seq library was prepared using a strand-specific method and you supplied the correct STRAND_SPECIFICITY parameter when running CollectRnaSeqMatrics.

If your library was prepared with an non-strand-specific protocol then Picard CollectRnaSeqMatrics can't know whether a read is incorrectly mapped. If your reads are not from a strand-specific library you should specify "STRAND_SPECIFICITY=NONE" when running CollectRnaSeqMatrics and the output will report "0" for both "CORRECT_STRAND_READS" and "INCORRECT_STRAND_READS". See the Picard Metrics Definitions here.

[batman]

Thanks kmcarr,

Assuming I correctly specify the STRAND_SPECIFICITY, then does picard verify the correctness of strands against the annotation provided in REF_FLAT option? My current thinking is that it is impossible to verify the correctness of strandedness of any region not in this annotation.

[kmcarr]

Thanks kmcarr,

Assuming I correctly specify the STRAND_SPECIFICITY, then does picard verify the correctness of strands against the annotation provided in REF_FLAT option? My current thinking is that it is impossible to verify the correctness of strandedness of any region not in this annotation.

Sorry, I don't know what it does with reads not overlapping an annotated gene. You could query the developers or run some tests yourself using simulated data where the count of correctly, incorrectly and non-annotated alignments is known to see what the output looks like.