从bam文件中提取行

2024-04-19 03:27:58 发布

您现在位置:Python中文网/ 问答频道 /正文

我有一个bam文件,如下所示:

samtools view pingpon.forward.bam | head
K00311:84:HYCNTBBXX:1:1123:2909:4215    0   LQNS02000001.1:55-552   214 28M *   0   0   TCTAGTTCAACTGTAAATCATCCTGCCC    AAFFFJJJJJJJJJJJJJJJJJJJJJJJ    AS:i:-6 XS:i:-6 XN:i:0  XM:i:1  XO:i:0  XG:i:0  NM:i:1  MD:Z:9T18   YT:Z:UU
K00311:84:HYCNTBBXX:1:1123:2909:4215    0   LQNS02000001.1:55-552   214 28M *   0   0   TCTAGTTCAACTGTAAATCATCCTGCCC    AAFFFJJJJJJJJJJJJJJJJJJJJJJJ    AS:i:-6 XS:i:-6 XN:i:0  XM:i:1  XO:i:0  XG:i:0  NM:i:1  MD:Z:9T18   YT:Z:UU
K00311:84:HYCNTBBXX:1:1123:2909:4215    0   LQNS02000001.1:55-552   214 28M *   0   0   TCTAGTTCAACTGTAAATCATCCTGCCC    AAFFFJJJJJJJJJJJJJJJJJJJJJJJ    AS:i:-6 XS:i:-6 XN:i:0  XM:i:1  XO:i:0  XG:i:0  NM:i:1  MD:Z:9T18   YT:Z:UU
K00311:84:HYCNTBBXX:1:1123:2909:4215    0   LQNS02000001.1:55-552   214 28M *   0   0   TCTAGTTCAACTGTAAATCATCCTGCCC    AAFFFJJJJJJJJJJJJJJJJJJJJJJJ    AS:i:-6 XS:i:-6 XN:i:0  XM:i:1  XO:i:0  XG:i:0  NM:i:1  MD:Z:9T18   YT:Z:UU
K00311:84:HYCNTBBXX:1:1123:2909:4215    0   LQNS02000001.1:55-552   214 28M *   0   0   TCTAGTTCAACTGTAAATCATCCTGCCC    AAFFFJJJJJJJJJJJJJJJJJJJJJJJ    AS:i:-6 XS:i:-6 XN:i:0  XM:i:1  XO:i:0  XG:i:0  NM:i:1  MD:Z:9T18   YT:Z:UU

我还有另一个文件,其中包含我感兴趣的ID,如下所示:

K00311:84:HYCNTBBXX:1:2223:15798:5692
K00311:84:HYCNTBBXX:2:2211:11414:30696
K00311:84:HYCNTBBXX:2:2223:28879:41581

理想情况下,我希望从bam文件中提取以IDs文件中的IDs开头的行。目前我正在使用我编写的代码,但它不起作用。任何帮助都将不胜感激!谢谢

import pysam
import re


forward = pysam.AlignmentFile('pingpon.forward.bam', "rb")
reverse = pysam.AlignmentFile('pingpon.reverse.bam', "rb")

ids = open("IDs_results_bed_reverse.txt", "w")

for line in reverse:
        if re.match("(.*)(I|i)ds(.*)", line):
            print(line)

Tags: asmdytxsxmnmbamuu
2条回答

第一个with语句是读取id文件并创建包含条目的字典。最后一个with将从另一个文件中读取条目,并将正确的条目放入字典中

import re

regex = re.compile('[A-Z0-9]+:\d+:[A-Z0-9]+:\d+:\d+:\d+:\d+')

with open('id.bam') as file:
    ids = {}
    for line in file:
        if regex.match(line):
            temp = line.replace('\n', '')
            ids[temp] = []

print(ids)

with open('list.bam') as file:
    for line in file:
        if regex.match(line):
            temp = line.replace('\n', '').split(' ')
            if temp[0] in ids:
                ids[temp[0]].append(line.replace('\n', ''))

print(ids)

https://www.biostars.org/p/165090/ 有人问了一个类似的问题,在这里得到了回答

相关问题 更多 >