Сценарий snakemake имеет доступ к stdin / stdout для обработки потока

Для рабочего процесса Snakemake мне нужно манипулировать тегами во многих файлах BAM, и я хотел бы обработать их, передав их через скрипт ( используя директиву Snakemake script:). Конкретный способ, которым я делаю это, с обработкой потока pysam.

infile = pysam.AlignmentFile("-", "rb")
outfile = pysam.AlignmentFile("-", "wb", template=infile)

for s in infile:
    (flowcell, lane) = s.query_name.split(':')[0:2]
    rg_id = ".".join([flowcell, lane])
    s.set_tag('RG',rg_id,'Z')
    outfile.write(s)

Этот скрипт хорошо работает в автономном режиме, но я не смог понять, как его интегрировать с помощью snakemake script директивы. Я предпочитаю этот способ минимизировать использование ввода-вывода и оперативной памяти.

Редактировать: прибегают к прямой загрузке, чтобы исправить тег RG.

# parameters passed from snakemake
bam_file = snakemake.input[0]
fixed_bam_file = snakemake.output[0]

bamfile  = pysam.AlignmentFile(bam_file, "rb")
fixed_bamfile = pysam.AlignmentFile(fixed_bam_file, "wb", template = bamfile)

for i, read in enumerate(bamfile.fetch()):
    (flowcell, lane) = read.query_name.split(':')[0:2]
    rg_id = ".".join([flowcell, lane])
    read.set_tag('RG',rg_id,'Z')
    fixed_bamfile.write(read)
    if not (i % 100000):
        print("Updated the read group for {} reads to {}".format(i, rg_id))

bamfile.close()
fixed_bamfile.close()

РЕДАКТИРОВАТЬ: Змеиный Мейк run: а также shell: директивы устанавливают workdir: каталог, в то время как script: Директива действует относительно каталога, в котором был запущен Snakefile (сохраняя все аккуратно и аккуратно). Отсюда проблема размещения потокового процессора под script:,

2 ответа

С помощью shell вместо script директива:

rule all:
    input:
        expand('{sample}_edited.bam'), sample=['a', 'b', 'c']


rule somename:
    input:
        '{sample}.bam'
    output:
        '{sample}_edited.bam'
    shell:
        '''
        cat {input} > python edit_bam.py > {output}
        '''

@Кришан, кажется, вы уже нашли решение, и если так, то, может быть, хорошо опубликовать его в качестве ответа.

Кроме того, вы можете использовать объект {workflow} чтобы получить каталог Snakefile и оттуда построить путь к вашему скрипту Python. Если ваша структура каталогов:

./
├── Snakefile
├── data
│   └── sample.bam
└── scripts
    └── edit_bam.py

Змеиный файл может выглядеть так:

rule all:
    input:
        'test.tmp',

rule one:
    input:
        'sample.bam',
    output:
        'test.tmp',
    shell:
        r"""
        cat {input} \
        | {workflow.basedir}/scripts/edit_bam.py >  {output}
        """

Выполнено с snakemake -d data ...

Кажется workflow объект не задокументирован, но проверьте эту ветку. Есть ли способ получить полный путь к Snake-файлу внутри Snake-файла?

Другие вопросы по тегам