如何在Python/Snakemake中精准解析BAM文件的唯一样本名

需求描述 #

我有一系列样本文件，文件名格式示例如下：

• SC892138_CTGAAGCT-ACTCTGAG_L001_001.star_rg_added.sorted.dmark.bam
• SC892138_unisample_L001_001.star_rg_added.sorted.dmark.bam
• SC892155_CTGAAGCT-ACTCTGAG_L001_001.star_rg_added.sorted.dmark.bam

目前代码仅能解析出SC######部分，需要修改解析逻辑，提取出SCxxxxxx_xxxxxx-xxxxxx或SCxxxxxx_unisample格式的完整样本标识。

现有代码 #

import os
import glob
import itertools
import pandas
from collections import defaultdict

workdir = os.environ['PWD']

## ---- The parser may have to be customized for each run ---- ##
def parse_sampleID(filename):
    return filename.split('/')[-1].split('_')[0]


fastqs = glob.glob('/P/A/T/H/S/*star_rg_added.sorted.dmark.bam')

d = defaultdict(list)
for key, value in itertools.groupby(fastqs, parse_sampleID):
    d[key] += list(value)

# Need to modify sampleIDs to bams not R1/R2 files
sampleIDs = d.keys()

修改方案 #

方法一：字符串分割法（简单直接） #

利用文件名中_L001这个固定分隔符，分割后取前半部分即可得到完整样本标识：

def parse_sampleID(filename):
    # 提取文件名（去掉路径）
    basename = filename.split('/')[-1]
    # 按_L001分割，取第一个元素就是完整样本标识
    return basename.split('_L001')[0]

方法二：正则表达式法（更健壮） #

如果文件名格式存在微小变动风险，用正则表达式匹配更可靠，确保精准提取SC开头的完整标识：

import re

def parse_sampleID(filename):
    basename = filename.split('/')[-1]
    # 匹配SC+数字+下划线+后续内容，直到_L001之前
    match_result = re.match(r'(SC\d+_.+?)_L001', basename)
    if match_result:
        return match_result.group(1)
    #  fallback：如果匹配失败，返回原逻辑的SC部分
    return basename.split('_')[0]

修改后的完整代码示例 #

以正则表达式法为例：

import os
import glob
import itertools
import pandas
import re
from collections import defaultdict

workdir = os.environ['PWD']

## ---- The parser may have to be customized for each run ---- ##
def parse_sampleID(filename):
    basename = filename.split('/')[-1]
    match_result = re.match(r'(SC\d+_.+?)_L001', basename)
    if match_result:
        return match_result.group(1)
    return basename.split('_')[0]


fastqs = glob.glob('/P/A/T/H/S/*star_rg_added.sorted.dmark.bam')

d = defaultdict(list)
for key, value in itertools.groupby(fastqs, parse_sampleID):
    d[key] += list(value)

sampleIDs = d.keys()

内容的提问来源于stack exchange，提问作者Genetics