Snakemake переменное количество файлов - PullRequest
Новая
       43

Snakemake переменное количество файлов

3 голосов
/ 21 марта 2020

Я нахожусь в ситуации, когда я хотел бы разбить мой рабочий процесс на переменное количество кусков, которые я не знаю заранее. Возможно, проще всего объяснить проблему конкретикой:

Кто-то передал мне демультиплексированные файлы FASTQ, используя bcl2fastq с опцией no-lane-splitting. Я хотел бы разделить эти файлы по полосам, сопоставить каждую полосу по отдельности, а затем, наконец, собрать все заново. Тем не менее, я не знаю заранее количество полос.

В идеале, я хотел бы получить решение, подобное этому, 

rule split_fastq_file: (...)  # results in N FASTQ files
rule map_fastq_file: (...)  # do this N times
rule merge_bam_files: (...)  # merge the N BAM files

, но я не уверен, что это возможно. Функция expand требует, чтобы я знал число дорожек, и не могу понять, как можно было бы использовать подстановочные знаки для этого.

Я должен сказать, что я довольно новичок в Snakemake, и что я, возможно, совершенно не понял, как работает Snakemake. Мне потребовалось некоторое время, чтобы привыкнуть думать о вещах «вверх ногами», сосредоточившись на выходных файлах и затем работая в обратном направлении.

1 Ответ

4 голосов
/ 22 марта 2020

Один из вариантов - использовать checkpoint при разбиении fastq, чтобы вы могли динамически пересмотреть DAG на более позднем этапе, чтобы получить результирующие дорожки.

Вот пошаговое руководство MWE:

  • Настройте и создайте пример файла fastq. 
# Requires Python 3.6+ for f-strings, Snakemake 5.4+ for checkpoints
import pathlib
import random

random.seed(1)

rule make_fastq:
    output:
        fastq = touch("input/{sample}.fastq")
  • Создайте случайное число дорожек от 1 до 9 со случайным идентификатором от 1 до 9. Обратите внимание, что мы объявляем это как checkpoint, а не как правило, чтобы впоследствии мы могли получить доступ к результату. Кроме того, мы объявляем выходные данные здесь как каталог, указывающий c для образца, чтобы мы могли позже использовать его, чтобы получить созданные дорожки. 
checkpoint split_fastq:
    input:
        fastq = rules.make_fastq.output.fastq
    output:
        lane_dir = directory("temp/split_fastq/{sample}")
    run:
        pathlib.Path(output.lane_dir).mkdir(exist_ok=True)
        n_lanes = random.randrange(1, 10)-
        lane_numbers = random.sample(range(1, 10), k = n_lanes)
        for lane_number in lane_numbers:
            path = pathlib.Path(output.lane_dir) / f"L00{lane_number}.fastq"
            path.touch()
  • Выполнить некоторую промежуточную обработку. 
rule map_fastq:
    input:
        fastq = "temp/split_fastq/{sample}/L00{lane_number}.fastq"
    output:
        bam = "temp/map_fastq/{sample}/L00{lane_number}.bam"
    run:
        bam = pathlib.Path(output.bam)
        bam.parent.mkdir(exist_ok=True)
        bam.touch()
  • Чтобы объединить все обработанные файлы, мы используем функцию ввода для доступа к дорожки, которые были созданы в split_fastq, так что мы можем сделать динамические c expand на них. Мы делаем expand на последнем правиле в цепочке промежуточных этапов обработки, в данном случае map_fastq, поэтому мы запрашиваем правильные входные данные. 
def get_bams(wildcards):
    lane_dir = checkpoints.split_fastq.get(**wildcards).output[0]
    lane_numbers = glob_wildcards(f"{lane_dir}/L00{{lane_number}}.fastq").lane_number
    bams = expand(rules.map_fastq.output.bam, **wildcards, lane_number=lane_numbers)
    return bams
  • Это теперь функция ввода дает нам легкий доступ к файлам bam, которые мы будем sh объединять, сколько бы их ни было, и как бы их ни называли. 
rule merge_bam:
    input:
        get_bams
    output:
        bam = "temp/merge_bam/{sample}.bam"
    shell:
        "cat {input} > {output.bam}"

Этот пример выполняется и с random.seed(1) создает три полосы (l001, l002 и l005).

Если вы не хотите использовать checkpoint, я думаю, вы могли бы добиться чего-то подобного, создав Функция ввода для merge_bam, которая открывает исходный ввод fastq, сканирует прочитанные имена для информации о дорожке и прогнозирует, какими должны быть входные файлы. Однако это кажется менее надежным.

...