💻

Easyfigでゲノムの比較結果を可視化（visualize synteny patterns using Easyfig）

に公開
Python
Bioinformatics
ゲノム解析
tech

はじめに

Genbank形式のアノテーションのついたファイルを使ってゲノムの比較結果を可視化してみます。
EasyfigをCLIで使用する例です。

動作確認 (Hardware)

-MacBook Pro
-チップ Apple M1 (Rosetta2 インストール済み)
-macOS Sonoma 14.7.5
-メモリ 16GB

動作確認 (Software)

-Python 3.12.10
-Easyfig.py Version 3.0.0
-blast 2.16.0

手順

1.データの用意

Nature Scientific Reportsに掲載の、こちらの文献を例にとります。
Comparative genomic analysis of Colistin resistant Escherichia coli isolated from pigs, a human and wastewater on colistin withdrawn pig farm | Scientific Reports

こちらの文献のFigure 3に相当する図を描いてみます。

1-1. Fastaファイルの入手

以前の記事に記した方法で、計6つのプラスミドのfastaファイルをダウンロードしました。

1-2. アノテーションを付与

アノテーションに使用したデータベースのバージョンなどの条件を各株で揃える目的で、自前でアノテーションを付与し直します。
今回は大腸菌（バクテリア）のプラスミドなので、baktaをローカル環境で動かしてアノテーションを付与しました。

baktaを使用すると拡張子が「.gbff」のGenbank形式のファイルが得られました。


2.プログラムの入手・インストール

2-1. Easyfigのプログラムファイル（Easyfig.py）の入手

Easyfigのメインサイトは以下のリンクです。

ただし、こちらの「DOWNLOADS」のメニューにあるEasyfig.pyのスクリプトは、Pythonのバージョン2対応版です。
Python2.7の仮想環境を作ってその中で動かすこともできますが、以下のGithubのリポジトリにPython3対応版のプログラムがあるので、こちらからEasyfig.pyをダウンロードして使用します。

コマンドを確認します。

shell
python path_to/Easyfig.py --help


2-2. blastのインストール

blastが入っていない場合はインストールしてPATHを通しておきます。
conda(mamba)を使用する場合の一例↓

shell
mamba install bioconda::blast #blastのインストール
which blastn #blastnにPATHが通っているか確認
which makeblastdb #makeblastdbにPATHが通っているか確認

3.描図

描図したいファイルを集めたディレクトリ上でEasyfigを実行していきます。

shell
cd path_to/genbank_files_dir

まずはシンプルに、株の並び順が文献と同じになるように描いてみます。
コマンドで指定したファイルの順に、図では上から表示されます。
アウトプットの画像ファイルの形式は、デフォルトではbmpファイルです。

shell
python path_to/Easyfig.py \
    -o test1.bmp \
    pCPA1200-IncI2.gbff \
    pCPWW7-IncI2.gbff \
    pHNSHP45-IncI2.gbff \
    pmcr1-IncI2.gbff \
    pCP52E-IncX4.gbff \
    pCSZ4-IncX4.gbff

test1



オプションのコマンドを追加します。

shell
python path_to/Easyfig.py \
   -o test2.bmp \
   -f1 T \ #blast検索の相同性のcolour gradientを図に表示
   -f2 5000 \ #scale barを5000bp（5Kbp）で図に表示
   -i 30 \ #blast検索のminimum identity valueを30%に設定
   -ann_height 100 \ #annotation（の図形）の高さを50ピクセル（デフォルト）から100ピクセルに変更
   -blast_height 200 \ #blast hitの高さを100ピクセル（デフォルト）から200ピクセルに変更
   -f CDS 255 165 0 arrow \ #CDS featureをOrange（RGB値：255 165 0）の矢印で表示
   pCPA1200-IncI2.gbff \
   pCPWW7-IncI2.gbff \
   pHNSHP45-IncI2.gbff \
   pmcr1-IncI2.gbff \
   pCP52E-IncX4.gbff \
   pCSZ4-IncX4.gbff

test2


ファイル名の後に整数値を2つ記述することで、配列中の描図する範囲を指定することができます。

描図する範囲を調べるには、Genbank形式のファイルから直接調べて決めてもよいですし、以前の記事に記した方法で、一旦CDS feature内のqualifierの情報を表形式で取得しておいてその表を見て範囲を絞るのに使う方法もあります。


各々の株の描図する範囲を指定して、さらに上下にgene qualifierの情報を表示させてみます。

shell
python path_to/Easyfig.py \
    -o test3.bmp \
    -f1 T \ #blast検索の相同性のcolour gradientを図に表示
    -f2 5000 \ #scale barを5000bp（5Kbp）で図に表示
    -i 30 \ #blast検索のminimum identity valueを30%に設定
    -ann_height 100 \ #annotation（の図形）の高さを50ピクセル（デフォルト）から100ピクセルに変更
    -blast_height 200 \ #blast hitの高さを100ピクセル（デフォルト）から200ピクセルに変更
    -f CDS 255 165 0 arrow \ #CDS featureをOrange（RGB値：255 165 0）の矢印で表示
    -legend both \ #qualifierの値を上下（図の一番上と一番下）の株で表示
    -leg_name gene \ #表示するqualifierをgeneに指定
    pCPA1200-IncI2.gbff 10996 22421 \
    pCPWW7-IncI2.gbff 11723 22696 \
    pHNSHP45-IncI2.gbff 15700 31805 \
    pmcr1-IncI2.gbff 13698 28000 \
    pCP52E-IncX4.gbff 19317 32301 \
    pCSZ4-IncX4.gbff 1424 17493

test3



特定のCDSに「-f」オプションで指定したものとは別の色を指定していきます。
個々のGenbank形式のファイルの特定のCDS featureに「colour」というqualifierを記述して色を指定します。

  • mcr-1に相当するCDSに Lime（RGB：0 255 0）
  • ISApl1に相当するCDSに Magenta（RGB：255 0 255）
  • pap2に相当するCDSに LightBlue（RGB：173 216 230）

図の一番上の株（pCPA1200-IncI2）のファイルへの記載例です。
mcr-1

ISApl1

pap2

これらの色の指定を各株のファイルに追記したのち再度描いてみます。

shell
python path_to/Easyfig.py \
    -o test4.bmp \
    -f1 T \ #blast検索の相同性のcolour gradientを図に表示
    -f2 5000 \ #scale barを5000bp（5Kbp）で図に表示
    -i 30 \ #blast検索のminimum identity valueを30%に設定
    -ann_height 100 \ #annotation（の図形）の高さを50ピクセル（デフォルト）から100ピクセルに変更
    -blast_height 200 \ #blast hitの高さを100ピクセル（デフォルト）から200ピクセルに変更
    -f CDS 255 165 0 arrow \ #CDS featureをOrange（RGB値：255 165 0）の矢印で表示
    -legend both \ #qualifierの値を上下（図の一番上と一番下）の株で表示
    -leg_name gene \ #表示するqualifierをgeneに指定
    pCPA1200-IncI2.gbff 10996 22421 \
    pCPWW7-IncI2.gbff 11747 22696 \
    pHNSHP45-IncI2.gbff 15700 31150 \
    pmcr1-IncI2.gbff 13698 28095 \
    pCP52E-IncX4.gbff 19317 32404 \
    pCSZ4-IncX4.gbff 2692 17493

test4

CDSの凡例の表示を変える手段としては、以下が考えられます。

  • Genbankファイル上で「leg_name」オプションに指定しているqualifier（今回の場合は「gene」）の記載を変更する（下図）
  • 「legend」オプションを使わずに画像を出力して、画像編集ソフト等で文字列を追加する

エラーが出たときの対処

実行の際に、以下のようなエラーが2つ出たので対処法を記しておきます。

1.「/bin/sh: makeblastdb: command not found」

Easyfig.pyのコードを確認すると

os.path.exists()で調べているPATH上にmakeblastdbがないことが原因と思われたので、シンボリックリンクを作成しました。

shell
#/usr/local/bin/ にmakeblastdbのシンボリックリンクを作成
ln -s path_to/makeblastdb /usr/local/bin/makeblastdb

2.「TypeError: 'float' object cannot be interpreted as an integer」

Easyfig.pyのファイルをテキストエディタなどで開き、プログラム中に4箇所ある
「range(y, y+ht+1)」の部分を、range関数内の値がint型になるように
「range(int(y), int(y+ht+1))」に置換して保存後、再度実行しました。

おわりに

元の文献の図と比較すると、使用したアノテーションのプログラムやデータベースの違いなどによる部分はありますが、概ね近い図を描くことができたのではないかと思います。

ほんの一例ですが、EasyfigのCLIでの使い方の参考になれば幸いです。

参考

カラーコードのRGB表記の値については以下を参考にさせていただきました。

Discussion