How can I extract base on specific loci from all reads overlapping that?
Hi, all. I'm using samtools to extract reads from BAM which I intrested. For example, if I want to extract all reads overlapping chr1_10000 and I know samtools view in.bam chr1:10000-10000 can do that. But how can I know which base on chr1_10000 in each reads? (especially when there has INDEL in reads)
• 4,171 views
•
link
2 answers
using sam2tsv : http://lindenb.github.io/jvarkit/Sam2Tsv.html
$ java -jar dist/sam2tsv.jar ~/scf7180000354095.bam -R ~/scf_7180000354095.fasta
#READ_NAME FLAG CHROM READ_POS BASE QUAL REF_POS REF OP
HISEQ:348:H2YWCBCXX:1:1115:19355:88769 163 scf7180000354095 0 A D 3 A M
HISEQ:348:H2YWCBCXX:1:1115:19355:88769 163 scf7180000354095 1 C B 4 C M
HISEQ:348:H2YWCBCXX:1:1115:19355:88769 163 scf7180000354095 2 A D 5 A M
HISEQ:348:H2YWCBCXX:1:1115:19355:88769 163 scf7180000354095 3 T D 6 T M
HISEQ:348:H2YWCBCXX:1:1115:19355:88769 163 scf7180000354095 4 T D 7 T M
HISEQ:348:H2YWCBCXX:1:1115:19355:88769 163 scf7180000354095 5 G E 8 G M
HISEQ:348:H2YWCBCXX:1:1115:19355:88769 163 scf7180000354095 6 A H 9 A M
HISEQ:348:H2YWCBCXX:1:1115:19355:88769 163 scf7180000354095 7 A H 10 A M
HISEQ:348:H2YWCBCXX:1:1115:19355:88769 163 scf7180000354095 8 T H 11 T M
• 0 views
•
link
For four bases: ATG and C, gives counts base wise:
sh /<path to="" igvtools="">/igvtools count -w 1 --bases <filename>.bam <filename>.wig hg19
• 0 views
•
link
Log in to answer this question.