How can I add NH:i:1 tag for every uniquely mapped read from BWA

[adrian]

Hi:

I have a BWA aligned BAM file. This BWA generated BAM file does not have NH:i:1 tag.

My aim is to create a BAM file where read is mapped uniquely and add NH:i:1 tag at the end of read.

I can filter BWA generated BAM file based on XT:A:U tag.


However, I am not sure how to add NH:i:1 at end of every alignment read.

Is there a way to do in pysam. I don't know exactly how to add NH:i:1 tag at end of read has XT:A:U tag. I appreciate any help.

Thank you.


from the following:

AMSAG1_299:7:2206:6444:162983 1187 1 249246 23 75M = 249327 156 CTGCGGGTGAAGTTCAAAGTACTATTCGGGTTTTTTTTTTAAAGTTTGGCTGGGTGCAGCGGCTCACGCCTGTAA @@@FFFFDDDH?FBHIGIIDGHGIJIJJIJJGHIIJJGI>BBBDCE@A=6;@@?,38CA587;;;>><??A8:>C X0:i:1 X1:i:1 XA:Z:16,-90270336,75M,1; MD:Z:75 RG:Z:128215 XG:i:0 AM:i:0 NM:i:0 SM:i:23 XM:i:0 XO:i:0 XT:A:U <<<--


TO:

AMSAG1_299:7:2206:6444:162983 1187 1 249246 23 75M = 249327 156 CTGCGGGTGAAGTTCAAAGTACTATTCGGGTTTTTTTTTTAAAGTTTGGCTGGGTGCAGCGGCTCACGCCTGTAA @@@FFFFDDDH?FBHIGIIDGHGIJIJJIJJGHIIJJGI>BBBDCE@A=6;@@?,38CA587;;;>><??A8:>C X0:i:1 X1:i:1 XA:Z:16,-90270336,75M,1; MD:Z:75 RG:Z:128215 XG:i:0 AM:i:0 NM:i:0 SM:i:23 XM:i:0 XO:i:0 XT:A:U NH:i:1 <<<---

[maubp]

Someone has asked the same question on the pysam user group:

https://groups.google.com/forum/#!topic/pysam-user-group/Hj3pNw4ohEM

If you are going to send duplicate questions to multiple places, please include the links to avoid wasting people's time. Thank you.

[adrian]

Thank you!

[maubp]

Presumably (since you asked on the pysam group) you can program in Python?

Do you want to add

NH:i:1

to every single read, or should there be a test of some kind? e.g. Only do this to mapped reads?

Essentially open the input BAM file, loop over all the reads - modify the

.tags

attribute to change the tags, and output the modified read. If you can find an example using pysam to filter a BAM file that would be a good basis to start from.

(You don't have to use pysam - see https://github.com/peterjc/picobio/b..._strip_tags.py for a related task which uses plain Python and the samtools command to handle SAM/BAM conversion - but for more complicated work pysam would be worth learning)