RNA seq strandedness (reverse-strand)

Question

[ZuelTech]

Hi. I have an RNA seq data that is reverse-stranded. The library type is ISR. I want to count reads mapping to the genomes. I am new to bioinformatics. Can you check if my code is correct?

Plus strand reads (Read2 mapped to forward strand) ?

PLUS_COUNT=$(samtools view -f 0x80 -F 0x10 "$BAM" "$SEG" | wc -l)

Minus strand reads (Read2 mapped to reverse strand) ?

MINUS_COUNT=$(samtools view -f 0x80 -f 0x10 "$BAM" "$SEG" | wc -l)

[GenoMax]

Use a standard tool like featureCounts as noted here --> RNA-seq analysis

Answer 1

[michael.ante]

It depends on the alignment method you used. If you aligned against the transcriptome, you can do it that way you described. In case you aligned against the genome, you need to take into account on which DNA-strand the transcript is on. For that you can use RSeQC's infer-experiment tool.

[ZuelTech]

I actually used STAR tool to align it to viral genomes.

[michael.ante]

You can use IGV to look at your alignment using your reference and its gene-annotation. You'll see how your reads orient themselves in regard to the genes.

Check also GenoMax' comment.

[ZuelTech]

My purpose is to count reads mapping to the genomes. Was my code above in doing so correct?

Answer 2

[Shohel]

Hi, you can use hisat2 for alignment then Sam tools and then featurecount for coun matrix... Hisat2 code is hisat2 -x genome_index -1 reads_R1.fastq -2 reads_R2.fastq -S output.sam --rna-strandness RF -p 8