Hi guys,
I’m working with some RNA seq data and get some problems.
I mapped the paired-end reads on tophat and got the result in sam/bam format. Now I’m using samtools to analyze the result.
I want to extract the reads that map to more than one place in the genome, and this is my command line:
Samtools view –h –f 0x100 in.bam > out.sam
There are no output alignmens in the out.sam except the head, which means there are no multi-mapped reads
However, I’ve run my own program in perl and find that there’re lots of reads whose IDs appear more than twice in the sam file, which means .these read mapped more than one place in the genome. (Two paired-end reads have the same ID, so most of the time, a certain ID appearing twice means it’s a unique mapping)
I want to know whether my command line is right? I mean whether 0x100 is the right parameter to extract multi-mapped reads.
I’m working with some RNA seq data and get some problems.
I mapped the paired-end reads on tophat and got the result in sam/bam format. Now I’m using samtools to analyze the result.
I want to extract the reads that map to more than one place in the genome, and this is my command line:
Samtools view –h –f 0x100 in.bam > out.sam
There are no output alignmens in the out.sam except the head, which means there are no multi-mapped reads
However, I’ve run my own program in perl and find that there’re lots of reads whose IDs appear more than twice in the sam file, which means .these read mapped more than one place in the genome. (Two paired-end reads have the same ID, so most of the time, a certain ID appearing twice means it’s a unique mapping)
I want to know whether my command line is right? I mean whether 0x100 is the right parameter to extract multi-mapped reads.
Comment