CommonWL for practical use 1
作成したCWL ツール・ワークフローを実際に動かしてみた（その1）
次世代シークエンサーDRY解析教本 改訂第2版（以下、DRY解析教本2）の目次が出版社のページから公開された。 その目次を見ての通り、Level1（準備編）、Level2（実践編）、Level3（応用編）の三編構成で、第1版と大きな構成に変化はない。 大きくは、Level2にメタゲノム解析、アッセンブル解析が加わって、Level3が総取っ替えされたといえば適切だろうか。 ちなみに、この本のtwitterのハッシュタグは #NGS_DAT で、この本に関わるtweetは、このハッシュタグを付けて是非。
しかしながら、大きな変更点が実はそれ以外にある。
Common Workflow Language (CWL)の追加だ。
Level2の一番最後(p331-338)に付録「CWLがあれば、DRY解析はもう怖くない」がある。
CWLについて日本語でコンパクトに解説してくれて、非常に有用である。
しかし、実はすごいのはここだけでない。
Level2で紹介した全てのワークフローが、CWLで書かれたWorkflowとしてこの記事のサポートサイトで公開されている。
つまり、ここにあるCWLをgit cloneしてきて自分にとって必要なように書き換えさえすれば、CWLを使ったNGSデータ解析が割と手間が少なくできる。
残念ながら、HISAT2-Stringtieによるワークフローは、DRY解析教本2には収録されなかった。 そこで、このワークフローをCWLで記述して実際に動かすことを目標として、2019年12月16日(月)に行われた23rd Workflow Meetupにて、自分の課題として取り組んだ。 Meetup当日のlive codingに関しては別に書かれている（「実際に、指先一つで立ち上げる CWL ツール・ワークフロー作成環境を試してみた」）ので、そちらを参照。 このエントリは、その後の顛末を。
まず、live codingで書いてもらったのはHISAT2のmappingの部分からなので、その前のindexを作成する必要がある。
もっともメジャーな生物だとすでにindexを直接ダウンロードできたりするが、いつもそうであるとは限らない。
FASTA形式のゲノム配列ファイルをゲットしてそれに対してindexを作成する。
幸い、pitagora-networkプロジェクトでHISAT2のCWLを作成し公開してくれているので、それをgit cloneしてありがたく利用させていただく。
#!/bin/bash
git clone https://github.com/pitagora-network/pitagora-cwl
cwltool pitagora-cwl/tools/hisat2/index/hisat2_index.cwl \
--index_basename org1 --reference_fasta organism1.fa
無事、indexができたので、先日のlive codingで作成していただいたrepoをgit cloneしてくる。
そして、それをつかってHISAT2のmappingを実行する。
#!/bin/bash
git clone https://github.com/tom-tan/cwl-live-coding-2019-12-16
cwltool cwl-live-coding-2019-12-16/cwl/hisat2-stringtie_wf_se.cwl \
--annotation organism1.gtf --fastq fibroblast.fq.gz \
--hisat2_idx_basedir . --hisat2_idx_basename org1 --nthreads 6
実行中に、Samtools view が並列化されていないことに気が付いたので、pitagora-cwl/tools/samtools/sam2bam/samtools_sam2bam.cwlを以下のように改変して並列化可能なように。
cwlVersion: v1.0
class: CommandLineTool
hints:
DockerRequirement:
dockerPull: quay.io/biocontainers/samtools:1.9--h46bd0b3_0
baseCommand: ["samtools", "view", "-b"]
arguments:
- prefix: -o
valueFrom: $(runtime.outdir)/$(inputs.output_filename)
inputs:
output_filename:
label: "Filename to write output to"
doc: "Filename to write output to"
type: string
default: out.bam
nthreads:
type: int
default: 1
inputBinding:
prefix: -@
input_sam:
label: "SAM format input file"
doc: "SAM format input file"
type: File
inputBinding:
position: 50
outputs:
bamfile:
type: File
outputBinding:
glob: $(inputs.output_filename)
$namespaces:
s: https://schema.org/
edam: http://edamontology.org/
s:license: https://spdx.org/licenses/Apache-2.0
s:codeRepository: https://github.com/pitagora-network/pitagora-cwl
s:author:
- class: s:Person
s:identifier: https://orcid.org/0000-0003-3777-5945
s:email: mailto:inutano@gmail.com
s:name: Tazro Ohta
$schemas:
- https://schema.org/docs/schema_org_rdfa.html
- http://edamontology.org/EDAM_1.18.owl
その上で、cwl-live-coding-2019-12-16/cwl/hisat2-stringtie_wf_se.cwlのsamtools_sam2bamの部分をローカルなCWLを参照するように、
samtools_sam2bam:
run: samtools_sam2bam.cwl
に書き換えて実行すると、sam2bam(samtools view)も並列で動くようになった。
これを使っていくつかのFASTQに対して実行し、複数のGTFを作成することができた。
そこで調子に乗って、stringtie mergeでそれらのGTFを束ねるプログラムもCWLで、ともう一つ作成してもらったCWLを動かそうとした。しかし、
#!/bin/bash
cwltool cwl-live-coding-2019-12-16/cwl/stringtie-merge.cwl \
--nthreads 6 --guide organism1.gtf --o_name merged.gtf \
--mergelist_txt mergelist.txt
とすると Error: cannot find file 'Brain.gtf' and /var/lib/cwl/stg4116466b-c5a7-4b66-8aac-76fe5a0ef6d1/mergelist.txt does not look like GFF! というエラーとなった。
どうも、テキストファイルの中に書いたファイルというのがうまくCWL内から見られていないようで。
今日のところはここで力尽きました（その2に続く）。