Snakemake按样本配对状态选择执行对应规则的问题求助

问题场景 #

我用Snakemake构建数据分析流程，输入包含样本配对信息的DataFrame，样本表内容如下：

sample,fastq_1,fastq_2,replicate,paired_end  
A1,A1_R1.fq.gz,A1_R2.fq.gz,1,True  
A2,A2.fq.gz,,1,False

流程中设置了样本通配符：

samples = pd.read_csv('config/samplesheet_valid.csv', sep=",").set_index("sample", drop=False)
wildcard_constraints:
    sample_id = "|".join(samples.index)

我写了依赖通配符的get_paired函数，想根据样本的paired_end字段选择执行对应规则：

def get_paired(wildcards):
    if samples.loc[wildcards.sample_id, "paired_end"]:
        return True
    else:
        return False

if get_paired:
    rule name_sort_bam_paired:
        input:
            bam_file = '{BASE_FOLDER}/bowtie2/filter/{sample_id}.flT.sorted.bam',
        output:
            sorted_bam = temp('{BASE_FOLDER}/bowtie2/filter/{sample_id}.clN.name.sorted.bam'),
            out_bam = temp('{BASE_FOLDER}/bowtie2/filter/{sample_id}.clN.bam'),
            out_sort_bam = '{BASE_FOLDER}/bowtie2/filter/{sample_id}.clN.sorted.bam'
        params:
            samtools = config['tools']['samtools'],
            bampe_rm_orphan = 'scripts/bampe_rm_orphan.py',
            python = config['tools']['python']
        threads: lambda wildcards: 10 if 10 < config['threads'] else config['threads']
        shell:
            """
            {params.samtools} sort -n --threads {threads} \
            -o {output.sorted_bam} \
            -T {wildcards.BASE_FOLDER}/bowtie2/filter/{wildcards.sample_id}.mLb.clN.name.sorted \
            {input.bam_file}

            {params.python} {params.bampe_rm_orphan} \
            {output.sorted_bam} {output.out_bam} --only_fr_pairs

            {params.samtools} sort --threads {threads} \
            -o {output.out_sort_bam} \
            -T {wildcards.BASE_FOLDER}/bowtie2/filter/{wildcards.sample_id}.clN.sorted \
            {output.out_bam}
            """
else:
    rule name_sort_bam_singled:
        input:
            bam_file = '{BASE_FOLDER}/bowtie2/filter/{sample_id}.flT.sorted.bam',
        output:
            out_bam = temp('{BASE_FOLDER}/bowtie2/filter/{sample_id}.clN.bam'),
            out_sort_bam = '{BASE_FOLDER}/bowtie2/filter/{sample_id}.clN.sorted.bam'
        params:
            samtools = config['tools']['samtools']
        threads: lambda wildcards: 10 if 10 < config['threads'] else config['threads']
        shell:
            """
            cp {input.bam_file} {output.out_bam}
            
            {params.samtools} sort --threads {threads} \
            -o {output.sorted_bam} \
            -T {wildcards.BASE_FOLDER}/bowtie2/filter/{wildcards.sample_id}.clN.sorted \
            {input.out_bam}
            """

但运行时所有样本都只执行name_sort_bam_paired规则，无法区分配对/单端样本。

问题原因 #

1. 函数对象被误判为布尔值：if get_paired:中，get_paired是一个函数对象，在Python里非空的函数对象会被判定为True，所以这个条件永远成立，导致name_sort_bam_paired规则被加载，else分支的规则根本不会被定义。
2. Snakemake规则定义是静态的：规则在流程初始化阶段就会全部加载完成，而get_paired需要具体的wildcards值才能返回对应样本的配对状态，初始化阶段无法获取这些动态值，所以这种全局if-else的方式根本无法实现动态规则选择。

解决方法 #

方法1：使用规则的condition参数（推荐） #

Snakemake支持给规则添加condition参数，该参数接受一个依赖wildcards的函数，返回布尔值，只有当返回True时，该规则才会被用于生成对应输出。

修改后的代码：

def get_paired(wildcards):
    return samples.loc[wildcards.sample_id, "paired_end"]

rule name_sort_bam_paired:
    condition: get_paired
    input:
        bam_file = '{BASE_FOLDER}/bowtie2/filter/{sample_id}.flT.sorted.bam',
    output:
        sorted_bam = temp('{BASE_FOLDER}/bowtie2/filter/{sample_id}.clN.name.sorted.bam'),
        out_bam = temp('{BASE_FOLDER}/bowtie2/filter/{sample_id}.clN.bam'),
        out_sort_bam = '{BASE_FOLDER}/bowtie2/filter/{sample_id}.clN.sorted.bam'
    params:
        samtools = config['tools']['samtools'],
        bampe_rm_orphan = 'scripts/bampe_rm_orphan.py',
        python = config['tools']['python']
    threads: lambda wildcards: 10 if 10 < config['threads'] else config['threads']
    shell:
        """
        {params.samtools} sort -n --threads {threads} \
        -o {output.sorted_bam} \
        -T {wildcards.BASE_FOLDER}/bowtie2/filter/{wildcards.sample_id}.mLb.clN.name.sorted \
        {input.bam_file}

        {params.python} {params.bampe_rm_orphan} \
        {output.sorted_bam} {output.out_bam} --only_fr_pairs

        {params.samtools} sort --threads {threads} \
        -o {output.out_sort_bam} \
        -T {wildcards.BASE_FOLDER}/bowtie2/filter/{wildcards.sample_id}.clN.sorted \
        {output.out_bam}
        """

rule name_sort_bam_singled:
    condition: lambda wildcards: not get_paired(wildcards)
    input:
        bam_file = '{BASE_FOLDER}/bowtie2/filter/{sample_id}.flT.sorted.bam',
    output:
        out_bam = temp('{BASE_FOLDER}/bowtie2/filter/{sample_id}.clN.bam'),
        out_sort_bam = '{BASE_FOLDER}/bowtie2/filter/{sample_id}.clN.sorted.bam'
    params:
        samtools = config['tools']['samtools']
    threads: lambda wildcards: 10 if 10 < config['threads'] else config['threads']
    shell:
        """
        cp {input.bam_file} {output.out_bam}
        
        {params.samtools} sort --threads {threads} \
        -o {output.out_sort_bam} \
        -T {wildcards.BASE_FOLDER}/bowtie2/filter/{wildcards.sample_id}.clN.sorted \
        {output.out_bam}
        """

注意：

• 两个规则的输出文件路径必须完全一致，这样Snakemake会根据condition判断用哪个规则生成该输出。
• 修复了原单端规则shell中的变量错误：原代码里output.sorted_bam不存在，改为output.out_sort_bam；input.out_bam改为output.out_bam。

方法2：使用统一规则+分支逻辑 #

如果不想拆分两个规则，可以在单个规则的shell脚本里根据样本配对状态执行不同逻辑：

def get_paired(wildcards):
    return samples.loc[wildcards.sample_id, "paired_end"]

rule name_sort_bam:
    input:
        bam_file = '{BASE_FOLDER}/bowtie2/filter/{sample_id}.flT.sorted.bam',
    output:
        sorted_bam = temp('{BASE_FOLDER}/bowtie2/filter/{sample_id}.clN.name.sorted.bam'),  # 仅配对样本使用
        out_bam = temp('{BASE_FOLDER}/bowtie2/filter/{sample_id}.clN.bam'),
        out_sort_bam = '{BASE_FOLDER}/bowtie2/filter/{sample_id}.clN.sorted.bam'
    params:
        samtools = config['tools']['samtools'],
        bampe_rm_orphan = 'scripts/bampe_rm_orphan.py',
        python = config['tools']['python'],
        is_paired = get_paired
    threads: lambda wildcards: 10 if 10 < config['threads'] else config['threads']
    shell:
        """
        if [ "{params.is_paired}" = "True" ]; then
            {params.samtools} sort -n --threads {threads} \
            -o {output.sorted_bam} \
            -T {wildcards.BASE_FOLDER}/bowtie2/filter/{wildcards.sample_id}.mLb.clN.name.sorted \
            {input.bam_file}

            {params.python} {params.bampe_rm_orphan} \
            {output.sorted_bam} {output.out_bam} --only_fr_pairs
        else
            cp {input.bam_file} {output.out_bam}
        fi

        {params.samtools} sort --threads {threads} \
        -o {output.out_sort_bam} \
        -T {wildcards.BASE_FOLDER}/bowtie2/filter/{wildcards.sample_id}.clN.sorted \
        {output.out_bam}
        """

这种方式把配对和单端的公共步骤（最后一步排序）抽出来，前面的分支根据is_paired参数判断执行逻辑。

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