python - 从 bam 文件中的 Cigar 字符串计算每个外显子的大小

Question

我有一个 bam 文件，我想计算每个外显子的大小：例如，我有 6 根雪茄串

1221M63I1145M
727M
288M
50M1I193M1I4M3D101M1I93M1D266M1D403M3I16M1D87M63I518M
383M26D1903M1I439M213N11M2D55M2I12M4D113M46N12M3D3M2D10M1D6M3I21M9D8M1D17M1D144M1D23M1I45M1I131M1I131M1I131M1D5M2D192M7I131M1D33M2I72M1D421M3D107M
6M3I21M9D8M1D17M1D144M1D23M1I45M1I131M1I131M1I131M1D5M2D192M7I135M2I3M3I27M2I72M1D421M3D112M1D22M5S

我想要输出：

exon: 2426  
exon: 727
exon: 288
exon: 1800
exon: 2726
exon: 193
exon: 1659
exon: 1667

这是我的脚本，但它只适用于前三个雪茄圈，我不知道如何才能得到正确的输出

import sys
import re
l=0
for line in sys.stdin:
    line = line.split("\t")
    t=re.split('([A-Z])',line[5])
    for i in range(0,len(t)-1,2) :
          if t[i+1] == 'N' :
             print("exon",l)
             l=0
          else:
              if t[i+1] == 'M' and len(t) <= 3:
                l=int(t[i])
              else:
                 if (t[i+1] == 'M' or t[i+1] == 'I') and len(t) > 3 :
                   l+=int(t[i])
    print("exon",l)

我的脚本的输出：

('exon', 2429)
('exon', 727)
('exon', 288)
('exon', 2088)
('exon', 4814)
('exon', 193)
('exon', 1659)
('exon', 3326)
('exon', 5195)
('exon', 85)

Answer 1

看起来您只是从每行中提取所有数字组并将它们相加。如果是这样的话：

import re
data = """1221M63I1145M
727M
288M
50M1I193M1I4M3D101M1I93M1D266M1D403M3I16M1D87M63I518M
383M26D1903M1I439M213N11M2D55M2I12M4D113M46N12M3D3M2D10M1D6M3I21M9D8M1D17M1D144M1D23M1I45M1I131M1I131M1I131M1D5M2D192M7I131M1D33M2I72M1D421M3D107M
6M3I21M9D8M1D17M1D144M1D23M1I45M1I131M1I131M1I131M1D5M2D192M7I135M2I3M3I27M2I72M1D421M3D112M1D22M5S"""

for line in data.split():
    numbers = re.findall(r'(\d+)', line)
    total = sum(int(n) for n in numbers)
    print(f'exon: {total}')

输出：

exon: 2429
exon: 727
exon: 288
exon: 1806
exon: 4895
exon: 1692