Hi
I have mapped bisulfite converted HiSeq reads separately to in silico converted forward & reverse strands so have files
I now would like a count of the reads which map to both and only both of the strands i.e. those which are in both forward & reverse bam files but not in only one of them.
Any ideas?
I was thinking of using either:
but can't work out how this would help
Could another option could be to use bedtools to convert both bam to bed & then intersect these? Something like:
Would this work?
Thanks
I have mapped bisulfite converted HiSeq reads separately to in silico converted forward & reverse strands so have files
Code:
sample.forward.bam sample.reverse.bam
Any ideas?
I was thinking of using either:
Code:
samtools merge out.bam sample.forward.bam sample.reverse.bam
Could another option could be to use bedtools to convert both bam to bed & then intersect these? Something like:
Code:
bamToBed -i sample.forward.bam > for.reads.bed bamToBed -i sample.reverse.bam > rev.reads.bed intersectBed -a .for.reads.bed -b rev.reads.bed -c
Thanks
Comment