Возможно ли иметь дополнительный выходной файл в Snakemake? - PullRequest
2 голосов
/ 30 мая 2019

Я пишу правило создания змеи, которое будет обрабатывать выполнение fastq-обрезки для данных секвенирования одного или парного конца. Если данные сопряжены, должно быть два выходных файла, если один конец, должен быть один.

Правило, которое я написал, работает до сих пор, однако у меня нет его, что вторая пара обрезанных является выходным файлом. Это означает, что snakemake не проверяет, существует ли этот файл. Он выдаст, но не проверяет, возможно ли иметь дополнительный вывод?

    input:
    #get the value in the fast1 column
        fastq_file = lambda wildcards: return_fastq(wildcards.fastq_name,wildcards.unit,first_pair = True)
    output:
        out_fastqc = config["fastp_trimmed_output_folder"] + "{unit}/{fastq_name}_trimmed.fastq.gz",
        fastpjson = config["fastp_trimmed_output_folder"] + "{unit}/{fastq_name}_fastp.json",
        fastphtml = config["fastp_trimmed_output_folder"] + "{unit}/{fastq_name}_fastp.html"
    params:
        fastp_parameters = return_parsed_extra_params(config['fastp_parameters']),
        fastq_file2 = lambda wildcards: return_fastq(wildcards.fastq_name,wildcards.unit,first_pair = False),
        out_fastqc2 = lambda wildcards: return_fastq2_name(wildcards.fastq_name,wildcards.unit),
        fastpjson = config["fastp_trimmed_output_folder"] + "{unit}/{fastq_name}_fastp.json",
        fastphtml = config["fastp_trimmed_output_folder"] + "{unit}/{fastq_name}_fastp.html"
    run:
        if config["end_type"] == "se":
            shell("{config[fastp_path]} -i {input.fastq_file} -o {output.out_fastqc} --json {output.fastpjson} --html {output.fastphtml} {params.fastp_parameters}")
        if config["end_type"] == "pe":
            shell("{config[fastp_path]} --in1 {input.fastq_file} --in2 {params.fastq_file2} --out1 {output.out_fastqc} --out2  {params.out_fastqc2} --json {output.fastpjson} --html {output.fastphtml} {params.fastp_parameters}")

Цель состоит в том, чтобы out_fastqc2 был включен в качестве необязательного вывода правила, чтобы snakemake проверил, существует ли он, и, если нет, соответствующим образом выдаст мне ошибку.

Если Snakemake не разрешает дополнительные выходы, я мог бы просто разделить на два правила, но это не совсем то, что я хотел бы.

Ответы [ 2 ]

0 голосов
/ 03 июня 2019

Поскольку кажется, что выбор между парным концом и односторонним в конфигурации является фиксированным, вы можете попытаться определить правило двумя способами в зависимости от конфигурации:

if config["end_type"] == "se":
    rule do_fastp:
        input:
        #get the value in the fast1 column
            fastq_file = lambda wildcards: return_fastq(wildcards.fastq_name, wildcards.unit, first_pair = True)
        output:
            out_fastqc = config["fastp_trimmed_output_folder"] + "{unit}/{fastq_name}_trimmed.fastq.gz",
            fastpjson = config["fastp_trimmed_output_folder"] + "{unit}/{fastq_name}_fastp.json",
            fastphtml = config["fastp_trimmed_output_folder"] + "{unit}/{fastq_name}_fastp.html"
        params:
            fastp_parameters = return_parsed_extra_params(config['fastp_parameters']),
        shell:
            """
            {config[fastp_path]} -i {input.fastq_file} -o {output.out_fastqc} \\
                --json {output.fastpjson} --html {output.fastphtml} \\
                {params.fastp_parameters}
            """
if config["end_type"] == "pe":
    rule do_fastp:
        input:
        #get the value in the fast1 column
            fastq_file1 = lambda wildcards: return_fastq(wildcards.fastq_name, wildcards.unit, first_pair = True),
            fastq_file2 = lambda wildcards: return_fastq(wildcards.fastq_name, wildcards.unit, first_pair = False)
        output:
            out_fastqc1 = config["fastp_trimmed_output_folder"] + "{unit}/{fastq_name}_trimmed_1.fastq.gz",
            out_fastqc2 = config["fastp_trimmed_output_folder"] + "{unit}/{fastq_name}_trimmed_2.fastq.gz",
            fastpjson = config["fastp_trimmed_output_folder"] + "{unit}/{fastq_name}_fastp.json",
            fastphtml = config["fastp_trimmed_output_folder"] + "{unit}/{fastq_name}_fastp.html"
        params:
            fastp_parameters = return_parsed_extra_params(config['fastp_parameters']),
        shell:
            """
            {config[fastp_path]} \\
                --in1 {input.fastq_file1} --in2 {params.fastq_file2} \\
                --out1 {output.out_fastqc1} --out2 {output.out_fastqc2} \\
                --json {output.fastpjson} --html {output.fastphtml} \\
                {params.fastp_parameters}
            """
0 голосов
/ 30 мая 2019

Посмотрите, как работает функция expand. Он вызывается на этапе, когда Snakemake создает DAG для зависимостей, и использует результат этой функции для создания списка файлов для раздела output.

Я бы предложил вам попробовать то же самое: составить список, который будет либо пустым, либо нет - зависит от условия.

Это решение будет работать, только если вы знаете, что вам нужно out_fastqc2 заранее (однако определение 2 правил с приоритетами делает то же самое). Если вы получаете информацию о необходимости out_fastqc2 только во время выполнения правила, это совершенно другой случай, когда вам нужны контрольные точки.

Ниже приведен код, иллюстрирующий мой подход: out_fastqc2 становится строкой, описывающей файл (если end_type настроен на "pe"), в противном случае становится пустым списком, который не меняет список выходы.

output:
    out_fastqc = config["fastp_trimmed_output_folder"] + "{unit}/{fastq_name}_trimmed.fastq.gz",
    out_fastqc2 = lambda wildcards: return_fastq2_name(wildcards.fastq_name,wildcards.unit) if config["end_type"] == "pe" else []
Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...