How to calculate GC content of reads that mapped to a specific gene?
1
1
Entering edit mode
9 months ago
bioinfo ▴ 150

Hello,

I aligned my data with STAR and I have the BAM files. I would like to calculate the GC content of reads that mapped to a specific gene. I have found this thread Gc Content From Bam but this does not seem to specify to do this just for a specific gene.

Thank you

STAR RNA-seq • 714 views
ADD COMMENT
4
Entering edit mode
9 months ago
bk11 ★ 2.7k

First create bam file for your gene-

samtools view -b input.bam "ChrX:10000-40000" > your_gene.bam

And then

samtools view  your_gene.bam | cut -f 10 | fold -w 1 | awk '($1=="G" || $1=="C") {N++;} END {print (N/(1.0*NR)*100);}'

OR you can do in using BBMap Suite

reformat.sh -Xmx2g in=your_gene.bam out=your_gene.fa | stats.sh -Xmx2g in=your_geneIN.fa
ADD COMMENT
0
Entering edit mode

Thank you so much. I used samtools and it worked perfectly.

ADD REPLY

Login before adding your answer.

Traffic: 1982 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6