BUSCOを使う

ゲノムのアセンブルに対する評価の一つに、必須遺伝子がどのぐらい見つかるかというものがあります。遺伝子が染色体に散在していると仮定して、あるはずの遺伝子が見つかる割合がゲノム全体に対するアセンブルされた長さの割合と一致するだろうという検証方法です。

これを実現するソフトウェアの一つがBUSCO (https://busco.ezlab.org/)です。そのバージョン5.3.1をdockerで使ってみました。

dockerイメージの取得

以下もそうですが、dockerを使うところはrootでないとだめです。

# docker search busco  # イメージを検索
INDEX       NAME                                      DESCRIPTION                                     STARS     OFFICIAL   AUTOMATED
docker.io   docker.io/chrishah/busco-docker           Docker image containing a full BUSCO install    1                    [OK]
docker.io   docker.io/ezlabgva/busco                  BUSCO - Benchmarking sets of Universal Sin...   1
docker.io   docker.io/robsyme/busco                   BUSCO - for  assessing genome assembly and...   1                    [OK]
(以下省略)
# docker pull docker.io/ezlabgva/busco:v5.3.1_cv1  # latestでないので:tagをつけないとだめ。使えるタグを調べる方法はないと思われる(イメージの作成者にあたるしかない)。

コンテナの調査

とりあえずコンテナ内の環境を調べるためにシェルを起動します。

# docker run -it --rm ezlabgva/busco:v5.3.1_cv1 bash
busco_user@3772446d2eee:/busco_wd$ id
uid=1000(busco_user) gid=1001(busco_user) groups=1001(busco_user),1000(busco_users)
busco_user@3772446d2eee:/busco_wd$ pwd
/busco_wd   # <- これが作業ディレクトリ
busco_user@3772446d2eee:/busco_wd$ busco -h  # 使い方を調べる
usage: busco -i [SEQUENCE_FILE] -l [LINEAGE] -o [OUTPUT_NAME] -m [MODE] [OTHER OPTIONS]
(以下省略)
busco_user@3772446d2eee:/busco_wd$ busco --list-datasets    # 使用できる遺伝子セット
(途中省略)
eukaryota_odb10   (真核生物)
     - viridiplantae_odb10  (緑色植物亜界)
         - embryophyta_odb10    (陸上植物)
             - eudicots_odb10    (真正双子葉類)
busco_user@3772446d2eee:/busco_wd$ exit   # --rmをつけているので自動的にコンテナが削除される。

BUSCOを実行

dockerの起動オプション

BUSCOの出力をローカルのファイルシステムに書き出すために、ホスト側のディレクトリをコンテナ内のディレクトリにマウントします(-v /home/takamasa/busco:/busco_wd)。
コンテナ内から解析対象のゲノム配列(/genome_dir/genome.fa)を見えるようにする(-v /genome_dir:/busco_wd/genome_dir)。
結果のファイルを一般ユーザーが自由にできるようにする(-u takamasa:takamasa)
時間がかかるのでバックグラウンドで実行する(-d)。パラメーターの設定がちゃんとできるまでは-dの代わりに-itとつけておくといいでしょう。

buscoの起動オプション

-f: 結果を上書きする。エラー終了のディレクトリをいちいち削除してから実行するのが面倒なので。
-i: 入力ファイル。コンテナ内でのパスであることに注意。
-o: 出力ディレクトリ。-hでは"Output folders and files will be labelled with this name."となっているが、実際にはこの名前のディレクトリが作成されて、その中にファイルがつくられる。
-m: 解析対象の配列の種類。RNA-Seqのアセンブルであればtranscriptomeを、アミノ酸配列であればproteinsを指定する。
-l: 解析対象の生物の系統。今回は双子葉植物を指定した。–list-databasesででてくる名前から適当なものを選ぶ

docker run --rm -d -u takamasa:takamasa \
    -v /home/takamasa/busco:/busco_wd \
    -v /genome_dir:/busco_wd/genome_dir \
    ezlabgva/busco:v5.3.1_cv1 \
    busco -f -i genome_dir/genome.fa -o out -m genome -l eudicots_odb10

380Mbpのゲノムで4.5時間かかりました。-cオプションで使用CPU数を指定できるようなので、これを増やすともう少し速くなるかもしれません。

結果の解釈

/home/takamasa/busco/out ディレクトリにあるshort_summary.specific.eudicots_odb10.out.txtを見る。

        ***** Results: *****

        C:56.2%[S:54.6%,D:1.6%],F:5.4%,M:38.4%,n:425
        239     Complete BUSCOs (C)
        232     Complete and single-copy BUSCOs (S)
        7       Complete and duplicated BUSCOs (D)
        23      Fragmented BUSCOs (F)
        163     Missing BUSCOs (M)
        425     Total BUSCO groups searched

調べた425個の遺伝子のうち、239個(232 + 7)がゲノムに見つかった。単純に言えば、アセンブルされた配列は全体の56%だよ、ということ。

run_eudicots_odb10

ここに詳細な結果がある。full_table.tsv と missing_busco_list.tsv は使用した遺伝子セットのうち解析対象に見つかったものと見つからなかったものが書かれている。
サブディレクトリのbusco_sequencesには見つかったアミノ酸配列がある。

ベンチマーク

ちなみに、シロイヌナズナのゲノムを調べてみると以下のようになりました。

        ***** Results: *****

        C:99.3%[S:98.6%,D:0.7%],F:0.0%,M:0.7%,n:425
        422     Complete BUSCOs (C)
        419     Complete and single-copy BUSCOs (S)
        3       Complete and duplicated BUSCOs (D)
        0       Fragmented BUSCOs (F)
        3       Missing BUSCOs (M)
        425     Total BUSCO groups searched