Snakemakeを使って再現性のある!スケーラブルな!可搬性のある解析! - ④

はじめに

※本記事はバイオインフォマティクス Advent Calendar 2019の20日目の記事です。

Snakemakeのインストールと概要
Snakemakeを使って再現性のある!スケーラブルな!可搬性のある解析! - ① - ばいばいバイオ

チュートリアルの実行
Snakemakeを使って再現性のある!スケーラブルな!可搬性のある解析! - ② - ばいばいバイオ

Snakemakeでの並列化
Snakemakeを使って再現性のある!スケーラブルな!可搬性のある解析! - ③ - ばいばいバイオ

①~③と、Advent Calendar に備え着々と準備を続けてきたのですが、
当日に全然時間がとれず今回は少し薄い内容となっております、どうかお許しを。

Snakemake vs Nextflow

バイオインフォ界隈のワークフローツールはSnakemakeだけでなく、CWL、Galaxy、Nextflowなんかが含まれてきます。 多種多様であれがいい、これがいい、そもそもワークフローはいらないだの色々と騒がれています。

少し古い引用ですが、NCBIのハッカソンでこんな議論が繰り広げられていたらしいです。

  • CWLは労働集約的すぎるため、ほとんどすべてのメンバーに却下された。 CWLの経験を持つ少数の人々は、扱いの困難さ、学習コストの高さを問題点として挙げていた。
  • SnakemakeはNextflowよりも柔軟性が低いという理由で却下された。多くのユーザーはPythonを使用していると考えていたが、そうではない人もいた。
  • Nextflowは、入力と出力を管理するだけの定義言語があらゆる言語のAPIで簡単にラップできるようになっていたため、採択された。
  • Jupyterは、まだすべての言語で扱えるものではない。将来的には魅力的な選択肢となるかもしれない。

当時はなんとも言えない理由でNextflowが採択されたようです。

参考:NCBI Hackathons discussions on Bioinformatics workflow engines

Snakemakeのいいところ

自分はCWL、Galaxyを少しかじった程度しかワークフローを触ったことはなく、pythonでごりごり書くことが多かったのですが、 巷ではやり始めたワークフローツールを使えばこれもっと楽にかけるんじゃね。。と思い、最近ではSnakemakeを信仰しています。

個人的に思うSnakemakeを使うべき人は下記の通りです。

  • Pythonが好き、Pythonでワークフローを記述したい
  • 小難しい文法を勉強したくない(学習コストをかけたくない)
  • プログラミング初学者が多いチームで使いたい
  • 数TBレベルの大規模処理は扱わない

ほかのワークフローツールと比べて、いろんな意味で軽量なツールかと思います。

GATK best practiceをSnakemakeで

dna-seq-gatk-variant-calling という便利なモノを見つけてしまいました。

これは、Snakemakeの作者であるJohannes Köster氏が作成したGATK best practiceのワークフローです。
彼はBiocondaプロジェクトのfounderや、Rust Bioの著者でもあるみたいです。

ワークフロー概略

Snakemakeのレポーティング機能は、下記のようなワークフロー図を出力してくれます。
少し画質が荒くて恐縮ですが、この図を張り付けておきます。

見てわかる通り、このワークフローはGATK best practiceに忠実なフローとなっています。
ワークフローの具体的な内容は本家を参考にしてください。

入力ファイル

units.tsvとして下記のように定義しています。
PLATFORMとか指定できるようになっていますが、レポーティング用途として使われているのみです。

sample unit platform fq1 fq2
A 1 ILLUMINA data/reads/a.chr21.1.fq data/reads/a.chr21.2.fq
B 1 ILLUMINA data/reads/b.chr21.1.fq data/reads/b.chr21.2.fq
B 2 ILLUMINA data/reads/b.chr21.1.fq

https://github.com/snakemake-workflows/dna-seq-gatk-variant-calling/blob/master/units.tsv

設定ファイル

config.yamlにて解析のオプションを設定します。
base recalibration を行うか否か、snv, indelごとのフィルタリング条件、trimmomaticやpicard, GATKのオプションが直感的に指定できるようになっています。

https://github.com/snakemake-workflows/dna-seq-gatk-variant-calling/blob/master/config.yaml

準備

このワークフローではBWAで使用するindex及びGATKが使用する.dictファイルが作成されていることが前提となっています。

それぞれリファレンスゲノムを指定して下記の通り実行することで作成できます。

bwa indexの実行

$ bwa index data/ref/genome.chr21.fa                                                                                         (snakemake) 

picard CreateSequenceDictionaryの実行

$ picard CreateSequenceDictionary R=data/ref/genome.chr21.fa  O=data/ref/genome.chr21.dict

実行

Snakemakeには便利なdry-run機能があります。
重い処理のワークフローを実行する前にはまずdry-runして概要が正しいか確認するようにしましょう。

$ snakemake --use-conda -n

実際にワークフローを実行する際は、--coresを指定して並列に処理を走らせることができます。

$ snakemake --use-conda --cores 12

中で使われている便利な機能たち

ログ出力

ログの出力はlogディレクティブに出力ファイル名を指定するだけ。
標準出力が吐き出されます。

明示的なパラメータを含む場合はshellディレクティブ内で指定してあげます。

rule abc:
    input: "input.txt"
    output: "output.txt"
    log: "logs/abc.log"
    shell: "somecommand --log {log} {input} {output}"

とっても簡単です。

参考:log-files

モジュール化

Snakemakeでは、includeディレクティブを使って処理のモジュール化を行うことができます。
ワークフローが大規模になってきた際には積極的に使っていきましょう。

使い方はいたって簡単。

下記のように含めたいSnakefileのパスをincludeディレクティブに指定するだけ。
モジュール内では通常のSnakefileと同じシンタックスが使用できます。

include: "path/to/other.snakefile"

参考:moduralization

wrapper

wrapperディレクティブを使うと、一般的な解析ツールであれば自分でコーディングせずにSnakemakeのワークフローの一部として使用できます。
例えばbwaだとこんな感じ。

rule bwa_mem:
  input:
      ref="data/genome.fa",
      sample=lambda wildcards: config["samples"][wildcards.sample]
  output:
      temp("mapped_reads/{sample}.bam")
  log:
      "logs/bwa_mem/{sample}.log"
  params:
      "-R '@RG\tID:{sample}\tSM:{sample}'"
  threads: 8
  wrapper:
      "0.15.3/bio/bwa/mem"

このディレクティブが指定されると、呼び出し時に要求されたバージョンのラッパーを自動的にダウンロードしてくれます。
--use-condaと組み合わせて使うと、wrapper内で必要なツールの仮想環境へのインストールも自動的に行ってくれます。
割と頻繁に更新されてるし簡単にオプションのカスタムできるし、とっても便利。

まだ数十ツールしかラインナップがないので優秀な皆さん、committerチャンスです。

参考: tool-wrappers

最後に

このワークフローはあくまでテンプレートとして使用して、カスタムして使うのが良いでしょう。
多言語対応という面ではNextflowなどにかなわないSnakemakeですが、Python愛好者には使いやすいことこの上なしです。

良いバイオインフォライフを。