EGAPxをローカルDockerで動かす

EGAPxをローカル環境で動かしたときの備忘録。

EGAPxはNCBIの真核生物ゲノムアノテーションパイプラインの外部実行版です。
ゲノム配列、taxid、RNA-seq、近縁種タンパク質などを使って、
GFFやタンパク質配列、CDS配列などを出してくれます。

NCBI RefSeq のアノテーションっぽい結果が自分の手元で出せて、
所感では他のソフトより断然いい感じのものができるが、かなり重めなソフトではあります。

大きめのゲノムで、しかもRNA-seqを入れると、普通に数日かかるし、途中でメモリ不足や時間切れに遭遇します。以下はローカルDockerで動かしたときのやり方と、詰まった時のTips。

公式ページ

GitHubはこちら https://github.com/ncbi/egapx

必要なもの

私はDocker executorで動かしました。

必要なものはだいたい以下。

  • Docker
  • Python3
  • Java 17以上
  • Nextflow
  • git
  • 十分なCPU、メモリ、ディスク（目安:32 thread以上、メモリ200 GB以上、ディスク100Gb以上）

EGAPxのデフォルト設定でも、huge jobは31 CPU / 200 GBになっている。
ただしこれは「重いtask 1個」の目安なので実際はもっと使うこともある。
ディスクもwork/以下に大量の中間ファイルができます。
3Gbのゲノムで中間ファイルが諸々100Gb以上できる。ひえ～。

インストール

まずEGAPxをclone。

git clone https://github.com/ncbi/egapx.git ~/egapx
cd ~/egapx

Python環境を作る。

python3 -m venv venv
source venv/bin/activate
pip install --upgrade pip
pip install -r requirements.txt

Nextflowはバージョン指定した方が無難。
私は23.10.1で動かしました。

export NXF_VER=23.10.1

Javaが古いとNextflowが動かないので確認。

java -version
nextflow -version
docker info

Dockerのpermissionで怒られる場合は、ユーザーがdocker groupに入っているか確認する。

id -nG

docker groupに入っているのに今のshellで反映されていない場合は、

newgrp docker

または、

sg docker -c 'コマンド'

で実行する。

入力yaml

EGAPxはyamlファイルに入力を書く。

最低限はゲノムFASTAとtaxid。
RNA-seqを入れる場合はshort_readsやlong_readsに書く。

例。

genome: /path/to/genome.fna
taxid: 12345
short_reads:
  - SRR0000001
  - SRR0000002
YAML

SRA accessionを指定するとEGAPx側で取得してくれる。
もちろんローカルFASTQを指定することもできる。

以下のようにNCBIの検索の時のクエリをそのまま入れても勝手に取得してくれる。

short_reads: txid12345[Organism] AND 50:350[ReadLength] AND (illumina[Platform] OR bgiseq[Platform]) AND biomol_rna[Properties]
YAML

タンパク配列のヒントファイルはtaxidから最適なセットを自動で取得されて使われる（もちろん指定もできる）。

Noncoding RNAのアノテーションも欲しい場合は、以下を追記。

cmsearch:
  enabled: true
trnascan:
  enabled: true
YAML

実行

基本形はこんな感じ。

python3 ~/egapx/ui/egapx.py egapx.yaml -e docker -o egapx.out -v

ただし、実際には時間がかかるのでnohupで投げる。

nohup python3 ~/egapx/ui/egapx.py egapx.yaml -e docker -o egapx.out -v > RunEGAP.nohup.out 2>&1 &

出力

完走するとegapx.out/以下に最終出力ができます。

よく使うのはこのへん。

egapx.out/complete.genomic.gff
egapx.out/complete.genomic.gtf
egapx.out/complete.genomic.fna
egapx.out/complete.proteins.faa
egapx.out/complete.cds.fna
egapx.out/complete.transcripts.fna
egapx.out/annotated_genome.asn

annotated_genome.asnはGenBank submissionで使えるらしい。

Tips

以下、実際に動かしていて効いたもの。

EGAPxの内部設定を確認したいときは、cloneしたrepoの中を見るとよい。
特にリソース設定はここ。

~/egapx/ui/assets/config/process_resources.config

ここにデフォルトのCPU、メモリ、time limit、retry回数が書かれている。

最初はRNA-seqなしで動かしてみる

まずはRNA-seqなしで動かしてみるのがおすすめ。

short_reads: []
long_reads: []
YAML

RNA-seqのヒントなしでも、かなり良いアノテーションができる。
RNA-seqは入れた方が常に良い、というものでもないらしい。
mappingが重いだけでなく、RNA evidenceが多すぎてGnomon/chainerあたりがめちゃくちゃ重くなることもある。

なので、最初からRNA-seq全部盛りで走らせるより、まずはミニマルな入力で走らせる。
それで質が悪くなさそうなら、その結果を使うのも全然あり。

関連ファイルのDLを使い回す

EGAPxはsupport dataやSRA readをNCBIのサーバーから自動でダウンロードする。
何度もランする場合、毎度DLし直したくないので、先にダウンロードしてそれを参照させる。

-lcでlocal cacheを指定できる。

python3 ~/egapx/ui/egapx.py egapx.yaml -dl -lc local_cache -e docker
python3 ~/egapx/ui/egapx.py egapx.yaml -lc local_cache -e docker -o egapx.out -v

resume

EGAPxはNextflowで動くので、途中で落ちてもresumeできる。

EGAPx実行後、だいたい以下のようなscriptができます。

egapx.out/nextflow/resume.sh

これを実行すると途中から再開できる。

bash egapx.out/nextflow/resume.sh

ただし、設定を変えて再開したい場合は注意。
yamlだけ直しても、すでに作られたrun_params.yaml側が使われていることがあります。

よく見るところ。

egapx.out/nextflow/run_params.yaml
egapx.out/nextflow/nextflow.log
egapx.out/nextflow/trace.txt

run_params.yamlを見て、今のresumeで本当に使われる設定を確認した方が良い。
resume時に使われるNextflow configは、resume.sh内のnextflow -C ...の行を見ると分かる。

STAR indexでメモリ不足になる

RNA-seqを入れるとSTARのindex作成が走る。
大きめのゲノムではここでメモリ不足になる。

STAR indexの工程はこのあたりに書かれている。

~/egapx/nf/subworkflows/ncbi/rnaseq_short/star_index/main.nf

このprocessはbig_jobなので、デフォルトではprocess_resources.configのbig_job設定が使われる。
エラーに必要メモリっぽい数字が出ることがあるので、それより大きい値を指定する。

yamlに追記する。

tasks:
  star_index:
    STAR: "--limitGenomeGenerateRAM 350000000000"
YAML

Nextflow側のメモリ制限も上げる必要がある。
追加configを作っておく。

process {
    withName: 'egapx:rnaseq_short_plane:star_index:build_index' {
        memory = 400.GB
        cpus = 8
        time = 12.h
    }
}
Bash

resume scriptのnextflow -Cにこのconfigを足して再開する。

miniprotでメモリ不足になる

タンパク質アライメントのところでminiprotがメモリ不足になることがある。

miniprotの工程はこのあたり。

~/egapx/nf/subworkflows/ncbi/target_proteins/miniprot/main.nf

run_miniprotにはhuge_jobとlong_jobが付いている。

エラー例。

[morecore] insufficient memory

単純にメモリを増やせば解決することもあるが、増やしてもダメな場合もある。
その場合は、タンパク質セットを減らす、chunkを小さくする、同時実行数を減らす、などを考える。

例。

tasks:
  miniprot:
    miniprot: "-t 8 -p 0.4 --outs=0.4"
  split_proteins:
    split_proteins: "-n 5000"
YAML

Nextflow側も調整する。

process {
    withName: 'egapx:target_proteins_plane:miniprot:run_miniprot' {
        memory = 800.GB
        cpus = 8
        time = 24.h
        maxForks = 1
    }
}
Bash

maxForks = 1にすると並列数は落ちるが、同時にメモリを食い尽くす事故は減る。

ProSplignはおすすめしない

EGAPxでは、protein alignmentにminiprotではなくProSplignを使う選択肢もあるが重いのであんまりおすすめしない。
（公式GitHubには、RNA-seqのヒントが少ない場合や、近縁種のprotein setが遠い場合に役立つ可能性があるとは書かれているが。）

run_prosplign_wnodeにもhuge_jobとlong_jobが付いていて、実際めちゃんこ時間がかかる（2周間で終わらなかった）。
まずはデフォルトのminiprotで動かす。
miniprotで結果が出たあと、どうしても必要ならProSplignを検討する、くらいでよいと思う。

なかなか終わらないときはchainerの16時間制限のせいかも

Gnomonのchainer stepが長時間走って、Nextflowのtime limitに当たることがある。

chainerは、RNA-seqやタンパク質のアライメントをもとに、ゲノム上でどのexon/intronをつないで遺伝子モデルにするかを探す工程。
STARやminiprotでできたヒントをただそのまま貼るのではなく、同じscaffold上の大量のevidenceを見ながら、Gnomonに渡す候補を組み立てている。

なので、ヒントが多い場所では重くなる。
特にRNA-seqのヒントが大量にあるscaffoldや、タンパク質アライメントが密集している領域では、組み合わせが増えてすごく時間がかかることがある。

EGAPxのデフォルト設定では、long_jobの制限時間は16時間。
chainerはこのlong_jobに入っている。

これは以下のファイルを見ると確認できる。

~/egapx/ui/assets/config/process_resources.config
~/egapx/nf/subworkflows/ncbi/gnomon/chainer_wnode/main.nf
process {
    errorStrategy = 'retry'
    maxRetries = 3

    withLabel: 'long_job' {
        time = 16.h
    }
}

つまり、16時間で終わらなかったrun_chainerはNextflowに止められて、自動でretryに入る。
maxRetries = 3なので、初回実行 + retry 3回で、最大4回試す。
同じchainer chunkがそれでも終わらないと、そこでworkflow全体がエラーになる。
つまりこのエラーで止まるまでに64h以上かかる（なんてこった）。
でも一回目で16h以内に終わらないなら、何度やり直しても大抵16hで終わらないのでここは見直すべき。

なかなか終わらないなと思ったら、まずnextflow.logを見る。

rg -n "run_chainer|Process exceeded running time limit|Execution is retried" egapx.out/nextflow/nextflow.log

16時間制限に当たっている場合は、こういう行が出る。

Process exceeded running time limit (16h)
NOTE: Process `...:chainer:run_chainer (4)` failed -- Execution is retried (1)

run_chainer (4)のような番号は、分割されたchainer jobの番号。
どのwork directoryで動いていたかもnextflow.logに出る。

workDir: /path/to/work/8b/9595f5...

そのwork directoryの中を見ると、実際の実行ログがある。

less /path/to/work/8b/9595f5.../.command.log

.command.logのタイムスタンプが更新され続けていて、request-startやrequest-stopが進んでいるなら、遅いがまだ進んでいる可能性が高い。
逆に、16時間付近でnextflow.logにProcess exceeded running time limit (16h)が出て、新しいwork directoryで同じrun_chainer番号が再 submitted されていたら、timeout retryに入っている。

解決策1 時間制限を伸ばす

例えばlong_jobの制限を伸ばすconfigを追加して無事終わることもある。

process {
    withLabel: 'long_job' {
        time = 48.h
    }
}

解決策2 ヒントを減らす

そもそもヒントが多すぎて時間のかかっているんだからヒントを減らせばいい。
RNAを使うのをやめれるならまずそれを試す。

RNA-seqを減らす以外に、タンパク質ヒントを減らすのも効くことがある。
EGAPxはtaxidから近そうなtarget protein setを自動で取ってくるが、このときに使う生物数をyamlファイルから絞れる。

proteins_best_n_orgs: 3
proteins_deny_taxids:
  - 1235
  - 4123

Docker内からNCBIに繋がらない

SRA取得などで、Docker container内からNCBIに繋がらないことがある。
ホストからは繋がるのにcontainer内から名前解決できない、というパターン。

その場合はDockerをhost networkで動かすと解決することがある。

Nextflow configの例。

docker.runOptions = "--network host --user $(id -u):$(id -g)"

--userも入れておくと、出力ファイルがroot所有になる事故を避けやすい。
EGAPxのDocker executorの元設定はここにある。

~/egapx/ui/assets/config/executor/docker.config

ここに--user $(id -u):$(id -g)が書かれている。
--network hostはそこに追加で足したもの。

コメント