如何用Snakemake按列表顺序配对ATAC与RNA序列生成libraries.csv

我需要为ATAC-seq和RNA-seq的SRR编号按顺序配对生成libraries.csv文件：现有2个ATAC-seq SRR编号和2个RNA-seq SRR编号，要求第一个ATAC SRR对应第一个RNA SRR，第二个对应第二个，最终生成2个{atac_srr}_{rna_srr}_libraries.csv文件。但当前Snakemake代码会生成所有可能的组合文件，请问如何修改实现顺序配对？

原代码 #

#Make 1 library.csv for a pair of ATAC and RNA seq accession numbers
def read_rna_accession():
    with open('rna_accessions.txt') as f:
        samples = [sample for sample in f.read().split('\n') if len(sample) > 0]  # Remove empty lines
        return samples

def read_atac_accession():
    with open('atac_accessions.txt') as f:
        samples = [sample for sample in f.read().split('\n') if len(sample) > 0]  # Remove empty lines
        return samples

# Read ATAC and RNA accession IDs
atac_SRRs = read_atac_accession()
rna_SRRs = read_rna_accession()


# Define all rule for generating libraries.csv
rule all:
    input:
        expand("{atac_srr}_{rna_srr}_libraries.csv", atac_srr=atac_SRRs, rna_srr=rna_SRRs)

# Rule to create libraries.csv
rule create_libraries_csv:
    output:
        "{atac_srr}_{rna_srr}_libraries.csv"
    run:
        atac_srr = wildcards.atac_srr
        rna_srr = wildcards.rna_srr
        with open(output[0], "w") as f:
            f.write("fastqs,sample,library_type\n")
            f.write(f"atac_seq/{atac_srr},{atac_srr},Chromatin Accessibility\n")
            f.write(f"rna_seq/{rna_srr},{rna_srr},Gene Expression\n")

修改方案 #

核心问题是原代码用expand生成了ATAC和RNA SRR的笛卡尔积（所有可能组合），我们需要改成按索引一一配对。修改后的代码如下：

# 为一对ATAC和RNA测序编号生成libraries.csv文件
def read_rna_accession():
    with open('rna_accessions.txt') as f:
        samples = [sample for sample in f.read().split('\n') if len(sample) > 0]  # 移除空行
        return samples

def read_atac_accession():
    with open('atac_accessions.txt') as f:
        samples = [sample for sample in f.read().split('\n') if len(sample) > 0]  # 移除空行
        return samples

# 读取ATAC和RNA的SRR编号
atac_SRRs = read_atac_accession()
rna_SRRs = read_rna_accession()

# 生成顺序配对的目标文件名列表（第一个ATAC对应第一个RNA，以此类推）
paired_files = [f"{atac_srr}_{rna_srr}_libraries.csv" for atac_srr, rna_srr in zip(atac_SRRs, rna_SRRs)]

# 定义所有需要生成的目标文件
rule all:
    input:
        paired_files

# 生成libraries.csv的规则
rule create_libraries_csv:
    output:
        "{atac_srr}_{rna_srr}_libraries.csv"
    run:
        atac_srr = wildcards.atac_srr
        rna_srr = wildcards.rna_srr
        with open(output[0], "w") as f:
            f.write("fastqs,sample,library_type\n")
            f.write(f"atac_seq/{atac_srr},{atac_srr},Chromatin Accessibility\n")
            f.write(f"rna_seq/{rna_srr},{rna_srr},Gene Expression\n")

关键改动说明 #

1. 替换笛卡尔积生成逻辑：用zip(atac_SRRs, rna_SRRs)将两个列表按索引一一配对，生成仅包含顺序配对的文件名列表paired_files，替代原代码中expand生成的全组合。
2. 保留核心生成规则：create_libraries_csv规则无需修改，因为它通过通配符依然能正确识别配对后的文件名中的ATAC和RNA SRR编号。
3. 注意事项：请确保atac_accessions.txt和rna_accessions.txt中的编号数量一致，zip会以较短列表的长度为准，自动截断多余的编号。

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