Thank you. However, I want to extract reads from Chr 10 with the given coordinates above.
Hello all,
I am extracting the reads from bam file which are falling in particular region using the following command.
samtools view -h input.bam "Chr10:18000-45500" > output.bam
However, now I want to extract the reads in sam manner but my reference is human genome primary assembly where the name of fasta reads starts like:
NC_000001.10 Homo sapiens chromosome 1, GRCh37.p13 Primary Assembly
How do I now modify command to do the sam task?
Thank you in advance.
1 answer
the fasta should be indexed ignoring everything after the first whitespace. So your command should be something like:
samtools view -h input.bam "NC_000001.10:123-456" > output.sam
try looking for some thing like: NC_000010 in header:)
$ grep ">" <input.fasta> (if your fasta is gizpped, use zgrep instead of grep). You would get all the headers in reference file. Then look for an entry that starts with NC_000010. Use that to extract the information. Command would look like:
samtools view -h input.bam "NC_000010.xx:123-456" > output.sam
Log in to answer this question.