snakemake, как построить цикл для двух независимых параметров - PullRequest
1 голос
/ 07 июня 2019

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

В случае, если для этого случая уже существует угроза, я был бы раднамек.Но до сих пор я не уверен, что правильные условия, чтобы искать то, что я хочу сделать.

Давайте предположим, что мой конвейер состоит из трех шагов.У меня есть набор образцов, которые я обрабатываю на каждом из этих трех этапов.Во втором шаге я добавляю дополнительный параметр к каждому образцу.На третьем шаге теперь я должен перебрать образцы и связанный параметр.Из-за этой структуры, я думаю, что это невозможно решить с помощью словарной структуры.

Чтобы представить пример, я сделал это упрощение моих файлов и правил: Файл конфигурации:

config.yaml

samples:
- a
- b
- c
- d
threshhold:
- 0.5
- 0.1
- 0.2

Схема моего змеиного файла.Он изображает точную структуру и наименование операций змеиного мейкера.По крайней мере, наименование упрощено.(Я добавил инструменты, которые я на самом деле использую в скобках, но я думаю, что это не обязательно для понимания.)

rule all: 
    input: 
      expand("{sample}.bam", sample=config["samples"]),
      expand("{sample}_{param}.bed", sample=config["samples"], param=config["threshhold"])

rule first: # (samtools view)
    input: 
       "{sample}.sam"
    output:
       "{sample}.bam"
    shell:
       "<somecommand> {input} {output}"

rule second: # ( macs2 callpeaks; Of course, there are multiple outputs but they only vary in their suffix))
    input:
       "{sample}.bam"
    output:
       "{sample}_{param}.bed"
    params:
       out_name="{sample}",
       threshhold="{param}"
    shell:
       "<somecommand> {input} -n {params.names} -q {params.threshhold}"

Итак, теперь у меня есть список файлов, выглядящих так:

  • a_0.5.bed
  • a_0.1.bed
  • a_0.2.bed
  • b_0.5.bed
  • b_0.1.bed
  • b_0.2.bed
  • ...

В моем третьем правиле я хочу сделать пересечения разных выборок с одним и тем же параметром.Например: a_0.5.bed x b_0.5.bed и c_0.5.bed x d_0.5.bed и получить вывод типа ab_0.5.bed, ab_0.1.bed , cd_0.5.bed ...

Моя первая попытка была такой:

rule all:
    input:
        expand("ab_{param}.bed", param=config["threshhold"])

rule intersect_2: # (bedtools intersect)
    input:
        a=expand("{sample_a}_{param}_peaks.narrowPeak", sample_a=config["samples"][0], param=config["threshhold"]),
        b=expand("{sample_b}_{param}_peaks.narrowPeak", sample_b=config["samples"][1], param=config["threshhold"])
    output:
        ab="intersect/ab_{param}.bed"
    params:
        threshhold="{param}"
    shell:
        "bedtools intersect -u -a {input.a} -b {input.b} > {output.ab}"

Ну, это не работает, потому что теперь входными данными являются все различные файлы параметров одновременно.

Я думаю, мне нужно больше и другая структура цикла здесь.Может быть, даже некоторые дополнительные петли Python вокруг правила или что-то?Но так как у меня совсем нет опыта программирования, и я только начинаю шаг за шагом разбираться с этими вещами, я уже не мог понять, с чего начать или какой цикл необходим для этого.

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

  • ab_0.5.bed
  • ba_0.5.bed
  • cb_0.5.bed
  • ca_0.5.bed
  • abc_0.5.bed
  • bca_0.5.bed
  • cba_0.5.bed

И эти комбинации тоже для всех остальных параметров.

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

РЕДАКТИРОВАТЬ: Может быть, поможет полностью реструктурированный файл конфигурации?Где образцы уже предварительно объединены?Может быть так: (предположим, что s1, s2 и т. Д. Означают реальное (и длинное) имя сэмпла)

config.yaml
samples_combinations:
- s1 : s2
- s3 : s2
- s3 : s1

Мне все еще нужно его переименовать ... Но я не оченьнравится идея.Моя цель состоит в том, чтобы создать что-то, что будет легко применимо и прямолинейно, без особого ручного уточнения, особенно потому, что в этом случае у вас может быть более трех сэмплов, которые нужно комбинировать несколькими способами.

1 Ответ

2 голосов
/ 07 июня 2019

Я думаю, ты почти понял. Вы можете достичь того, что вы хотите с простыми подстановочными знаками. Помните, что expand используется для создания списков. Сохраняйте это простым:

rule all:
    input:
        expand("intersect/ab_{param}.bed", param=config["threshhold"])

rule intersect_2: # (bedtools intersect)
    input:
        a="{sample1}_{param}_peaks.narrowPeak",
        b="{sample2}_{param}_peaks.narrowPeak"
    output:
        ab="intersect/{sample1}{sample2}_{param}.bed"
    shell:
        "bedtools intersect -u -a {input.a} -b {input.b} > {output.ab}"

Несколько замечаний по этому поводу:

  • Запутано использование подстановочных имен (или имен ввода / вывода), которые также являются реальными именами ваших семплов.
  • Довольно опасно ставить две или более подстановочные строки подряд без разделителей

Все вместе, я бы написал что-то вроде этого:

import itertools

sampleCombinations = [s[0]+"-"+s[1] for s in list(itertools.combinations(config["samples"],2))]

rule all:
    input:
        expand("intersect/{pairOfSamples}_{param}.bed", pairOfSamples=sampleCombinations, param=config["threshhold"])

rule intersect_2: # (bedtools intersect)
    input:
        s1="{sample1}_{param}_peaks.narrowPeak",
        s2="{sample2}_{param}_peaks.narrowPeak"
    output:
        "intersect/{sample1}-{sample2}_{param}.bed"
    shell:
        "bedtools intersect -u -a {input.s1} -b {input.s2} > {output}"

itertools.combinations создаст пары из ваших образцов, определенных в конфигурации.

Список sampleCombinations будет содержать все возможные пары выборок, разделенных дефисом (т. Е. «A-b», «b-c», «c-d» и т. Д.). Это может привести к поломке, если ваши имена семплов содержат дефисы, так как змеиный мастер не сможет хорошо восстановить символы подстановки ({sample1}-{sample2}). В этом случае выберите другой разделитель.

РЕДАКТИРОВАТЬ после комментариев ОП:

Извините за все мое правило, я забыл поставить вывод папки intersect.
Я также забыл _{param}_ номинал во входных файлах

Допустим, у вас есть файл конфигурации, подобный следующему:

samples:
- name: very_long_name_a
  shortName: a
- name: very_long_name_b
  shortName: b
- name: very_long_name_c
  shortName: c
threshhold:
- 0.5
- 0.1
- 0.2

тогда вы можете использовать что-то вроде этого:

import itertools

configfile: "config.yaml"

longNamesList = [s["name"] for s in config["samples"]]
shortNamesList = [s["shortName"] for s in config["samples"]]
shortToLong = {s["shortName"]:s["name"] for s in config["samples"]}

sampleCombinations = [s[0]+"-"+s[1] for s in list(itertools.combinations(shortNamesList,2))]

rule all:
        input: expand("intersect/{pairOfSamples}_{param}.bed", pairOfSamples=sampleCombinations, param=config["threshhold"])


rule intersect_2: # (bedtools intersect)
    input:
        s1=lambda wildcards: shortToLong[wildcards.sample1]+"_"+wildcards.param+"_peaks.narrowPeak",
        s2=lambda wildcards: shortToLong[wildcards.sample2]+"_"+wildcards.param+"_peaks.narrowPeak"
    output:
        "intersect/{sample1}-{sample2}_{param}.bed"
    shell:
        "bedtools intersect -u -a {input.s1} -b {input.s2} > {output}"
...