macでインフォマティクス

macでインフォマティクス

HTS (NGS) 関連のインフォマティクス情報についてまとめています。

特定の枝で正の選択が起きたか調べる aBSREL

2015 2018 2020 Molecular Biology and Evolution docker dispensable genes evolutionary model multiple sequence alignment (MSA) dN/dS phylogenetic tree viewer phylogenetic tree Inference

2026/01/14追記､ 1/23追記

 

　過去20年間にわたり、コドン置換モデルを用いた比較配列解析は、分子データから自然選択のシグネチャーを検出するための強力かつ一般的な手法へと磨き上げられてきた。多くの研究は、ω比で定量化された配列に対する選択圧をコドンサイト間および系統樹上の個々の枝間で変化させることを許容する「ブランチサイト」モデルの開発に焦点を当ててきた。本著者らは、このクラスの手法として、適応型ブランチサイトランダム効果尤度adaptive branch-site random effects likelihood（aBSREL）を開発し、提示する。この手法の重要な革新性は、情報理論的基準に基づいて選択される可変パラメトリック複雑度である。系統樹上の異なる枝に異なる複雑度のモデルを適用することで、aBSRELは、既存の最高クラスの手法に匹敵する、あるいはそれを上回る統計的性能を実現しながら、一桁高速に実行できる。シミュレーションデータ解析に基づき、多様化をもたらす正の選択の程度と強度をどの程度確実に検出できるかについてのガイドラインを示し、「ブランチサイト」モデルの最適なパラメトリック複雑度には自然な限界があることを示唆する。 8,893の真口動物の遺伝子アライメントを対象としたaBSREL解析では、典型的な遺伝子系統樹における分岐の80%以上が単一のω比モデルで適切にモデル化できることが示された。つまり、現在のモデルは不必要に複雑であるということである。しかしながら、モデル選択手順を用いてデータから同定される重要な分岐は比較的少数であり、進化の複雑さを正確にモデル化することが不可欠である。

 

HyPhy: Hypothesis testing using Phylogenies

http://hyphy.org

HyPhy Vision

http://vision.hyphy.org/

Datamonkey: HyPhyを使って配列データにおける進化シグネチャーを解析するためのオンライン解析プラットフォーム。インタラクティブ、API対応（最新のv3サーバーも利用可能）

Datamonkey Adaptive Evolution Server

aBSREL: adaptive Branch-Site Random Effects Likelihood

http://help.datamonkey.org/methods/absrel.html#example-workflow

 

aBSRELの使い方について簡単にみていきます。

インストール

ubuntu22とmacでテストした（macでも導入できたが動作せず）

HyPhy

HyPhy Vision

#conda
mamba create -n hyphy -y
conda activate hyphy
mamba install -c bioconda hyphy -y

#docker（未テスト）
git clone https://github.com/veg/hyphy.git
cd hyphy
docker build -t hyphy:latest .
docker run --rm -v [path-to-your-input-data]:/hyphy/data -it hyphy:latest

> hyphy absrel -h #-hなしで叩けばinterractiveモードで実行できる

 hyphy absrel

 

Analysis Description

--------------------

aBSREL (Adaptive branch-site random effects likelihood) uses an adaptive

random effects branch-site model framework to test whether each branch

has evolved under positive selection, using a procedure which infers an

optimal number of rate categories per branch. v2.2 adds support for

multiple-hit models. v2.3 adds support for SRV. v2.5 adds support for

ancestral state reconstruction, identification of sites contributing to

selection signal, and some diagnostics. v2.6 adds tree reports and ASCII

art

 

- __Requirements__: in-frame codon alignment and a phylogenetic tree

 

- __Citation__:  Less Is More: An Adaptive Branch-Site Random Effects Model for

Efficient Detection of Episodic Diversifying Selection (2015). Mol Biol

Evol 32 (5): 1342-1353.

 

- __Written by__: Sergei L Kosakovsky Pond, Ben Murrell, Steven Weaver and Temple iGEM /

UCSD viral evolution group

 

- __Contact Information__: spond@temple.edu

 

- __Analysis Version__: 2.6

 

[Advanced option] Bag of little bootstrap alignment resampling rate (permissible range = [0.01,1], default value = 1): blb: 1

 

>code => Universal

 

 

>Select a coding sequence alignment file (`/home/kazu/`) ^CError:

HyPhy killed by signal 2

 

Function call stack

1 :  [namespace = buWlJjnD] DataSet  ^dataset_name = ReadDataFile(PROMPT_FOR_FILE);

-------

2 :  [namespace = byBbVLgm] code_info=alignments.LoadGeneticCode(None);

-------

3 :  [namespace = zGIeuxOb] return alignments.ReadCodonDataSetFromPath(dataset_name,None);

-------

4 :  return alignments.LoadCodonDataFile(dataset_name,datafilter_name,alignments.ReadCodonDataSet(dataset_name));

-------

5 :  [namespace = absrel] codon_data_info=alignments.PromptForGeneticCodeAndAlignment(datasets,prefix+".codon_filter");

-------

6 :  [namespace = absrel] load_file({utility.getGlobalValue("terms.prefix"):"absrel",utility.getGlobalValue("terms.data.blb_subsample"):blb});

-------

7 :  namespace

 

Step 0.LoadFunctionLibrary("modules/shared-load-file.bf", /home/kazu/miniforge3/share/hyphy/TemplateBatchFiles/SelectionAnalyses/);

 

Step 1.load_file({utility.getGlobalValue("terms.prefix"):"absrel",utility.getGlobalValue("terms.data.blb_subsample"):blb});;

-------

 

Check errors.log for execution error details.

(base) kazu@DESKTOP-4D06TNK:~$ hyphy absrel -h

 

Available analysis command line options

---------------------------------------

Use --option VALUE syntax to invoke

If a [reqired] option is not provided on the command line, the analysis will prompt for its value

[conditionally required] options may or not be required based on the values of other options

 

blb

        [Advanced option] Bag of little bootstrap alignment resampling rate

        default value: 1.0 [computed at run time]

 

code

        Which genetic code should be used

        default value: Universal

 

alignment [required]

        An in-frame codon alignment in one of the formats supported by HyPhy

 

tree [conditionally required]

        A phylogenetic tree (optionally annotated with {})

        applies to: Please select a tree file for the data:

 

branches

        Branches to test

        default value: All

 

multiple-hits

        Include support for multiple nucleotide substitutions

        default value: None

 

srv

        Include synonymous rate variation

        default value: No

 

syn-rates

        The number alpha rate classes to include in the model [1-10, default 3]

        default value: absrel.synonymous_rate_classes [computed at run time]

 

output

        Write the resulting JSON to this file (default is to save to the same path as the alignment file + 'ABSREL.json')

        default value: absrel.codon_data_info[terms.json.json] [computed at run time]

 

intermediate-fits

        Use/save parameter estimates from 'initial-guess' model fits to a JSON file (default is not to save)

        default value: /dev/null

 

kill-zero-lengths

        Automatically delete internal zero-length branches for computational efficiency (will not affect results otherwise)

        default value: Yes

 

save-fit

        Save full adaptive aBSREL model fit to this file (default is not to save)

        default value: /dev/null

 

 

実行方法

in-frameのアラインメントデータが必要。翻訳領域のCDS配列を用意し、翻訳してアラインメントする。それから逆翻訳してコドンの位置を保ったままの整列結果を用意する（CDSのままアラインメントして使うとギャップが入りズレてしまう）。別途、系統推定結果も必要（同じ配列名で対応付ける）。

 

1、CDSの塩基配列を用意する。

 

2、翻訳。ここではembossのtranseqを使う。あるいはseqkitかgotransseqを使う。

#transeq (webサーバー)
transeq cds.fna protein.faa

#seqkit
seqkit translate cds.fna -o protein.faa

 

3、翻訳したAAを多重整列 - ここではmafftを使用

mafft --auto protein.faa > protein.aln

#精度を上げる
mafft --localpair --maxiterate 1000 protein.faa > protein.aln

#clustal omegaなら
clustalo -i protein.faa -o protein.aln --outfmt=fa --threads=4

 

（任意）MSAの確認 

#1 alvでターミナル上で可視化
#インストール
conda activate hyphy
mamba install -c bioconda alv
#実行
alv protein.aln

#2 あるいはmviewででターミナル上で可視化(-out clustalでclustal形式で可視化)
https://kazumaxneo.hatenablog.com/entry/2019/11/15/073000
#インストール
mamba install bioconda::mview -y
#ターミナル上で可視化
mview -in fasta -out clustal protein.aln
#HTML保存
mview -html head -coloring any -bold -css on -in fasta protein.aln > align.html

ALV

Mview

 

3、PAL2NAL（webサーバー）(github)を使って逆翻訳する。PAL2NALは忠実にコドンを塩基配列へ反映させるのでギャップもそのまま残る。

#condaで導入可能
conda activate hyphy
mamba install bioconda::pal2nal -y

#PAL2NALの実行
pal2nal.pl protein.aln cds.fna -output fasta > codon_align.fasta

> seqkit stats codon_align.fasta

file               format  type  num_seqs  sum_len  min_len  avg_len  max_len

codon_align.fasta  FASTA   DNA          4    3,960      990      990      990

(3の倍数であることも確認)

 

４、系統推定の実行。ここではIQtreeを使う。枝の長さが統計計算の基礎となるため、NJやMPではなく最尤（ML）法が推奨される。

#1 IQtreeの場合、latestは3.0.1
conda activate hyphy
mamba install bioconda::iqtree -y
iqtree3 -s codon_align.fasta -m MFP -B 1000 -T AUTO

#2 RAxML-NGの場合
conda activate hyphy
mamba install bioconda::raxml-ng -y
raxml-ng --msa codon_align.fasta --model GTR+G --prefix test --threads 2

#3 MEGA (CLI)の場合
conda activate hyphy
mamba install bioconda::mageck -y
#1 MEGA GUIでconfigファイル:maoファイルを作成。MLを設定、パラメータを決めて => save settigsでmaoファイル書きだし
#2 ML
megacc -a ml_template.mao -d codon_align.fasta -o my_ml_tree

 

5、HyPhyに実装されているaBSRELの実行。アラインメントとツリーファイルを指定する（http://help.datamonkey.org/methods/absrel.html）

conda activate hyphy
hyphy absrel --alignment codon_align.fasta --tree tree.nwk
  • --alignment   An in-frame codon alignment in one of the formats supported by HyPhy.
  • --tree   A phylogenetic tree (optionally annotated with {}).

=> .JSONファイルが書き出される。

 

 

6、HyPhy Visionにアクセスして分析、視覚化する: http://vision.hyphy.org/absrel

chromeでは動作しなかったのでsafariを使用した（M1 macbook air使用。os ventura）。

右上のLoadからJSONファイルを読み込ませる。

ω(Omega)は選択圧の強さ（dN/dS ratio）、１より大きいと正の選択の可能性。P値はその枝で正の選択が起こったという証拠が統計的に有意かどうか（一般的に0.05以下が有意）。4個の配列を含むこのテストデータでは有意な正の進化をしている枝はゼロだった（祖先の枝も含めて全て真っ青になっていて、5個のブランチが検定テスト対象となり、有意な枝はゼロ）。

 

ω分布の可視化: 各枝において、どの程度の割合のサイトがどのような選択圧を受けているか

初期読み込みデータ

 

 

6、iq-treeで推定祖先配列を取り出す｡例えばnode8で正の選択が起きていたら､node8の祖先配列を推定し､その配列を取り出してさらに調査することができる｡-asrですべての内部ノードの推定配列が出力される。

iqtree -s protein.aln -m LG+G4 -asr
=> 系統推論のファイルのほかに.stateファイルが書き出される

#正の進化の可能性が見られたnode8を抽出するなら､以下のようにする
echo ">Node8" > Node8.fasta
grep "^Node8:space:" nifH.aln.state | awk '{print $3}' | tr -d '\n' >> Node8.fasta

 

Datamonkeyでは、ガイドに従って適切なアプローチを選択することができる。

クリックして選択していくとaBSRELが適切と表示された。使いたいならaBSRELの部分をクリックする。

 

下の画面が表示される。ここでin-frameのアラインメントデータをアップすると、ABSREL側でそのまま系統推定も行われ、aBSRELによる計算結果付きの系統樹（上の図）が表示される。

 

追記

MEME: どのアミノ酸残基が変化したかについて分析できる｡

conda activate hyphy
hyphy meme --alignment codon_align.fasta --tree tree.nwk

meme.jsonが書き出される｡これをMEMEにアップロードする｡

http://vision.hyphy.org/MEME

MEME（Mixed Effects Model of Evolution）は、混合効果最大尤度法を用いて、個々のサイトがエピソード的な正の選択または多様化選択を受けているという仮説を検証する(HPより)。

 

BUSTED: 遺伝子全体として、正の選択圧がかかっている可能性を調査する｡

conda activate hyphy
hyphy busted --alignment codon_align.fasta --tree tree.nwk

BUSTED.jsonが書き出される｡これをBUSTEDにアップロードする｡

http://vision.hyphy.org/BUSTED

遺伝子が少なくとも1つの分岐において少なくとも1つの部位で正の選択を受けているか調べるため、遺伝子全体（部位特異的ではない）における正の選択の検定を行う(HPより)。

 

FEL: 純化選択（負の選択）」がかかっているサイトを調べる｡

conda activate hyphy
hyphy fel --alignment codon_align.fasta --tree tree.nwk

FEL.jsonが書き出される｡これをFELにアップロードする｡

http://vision.hyphy.org/FEL

FEL（Fixed EffectsLikelihood）は、最大尤度（ML）アプローチを用いて、与えられたコーディングアラインメントと対応する系統樹において、部位ごとに非同義置換率（dN）および同義置換率（dS）を推定する。この手法では、各部位に対する選択圧は系統樹全体にわたって一定であると仮定する(HPより)。

 

Datamonkeyのhelpより

  • アラインメントがフレーム内に収まっていること、つまり、早期STOPコドン（ミスアラインメントや機能しないコード配列などによるフレームシフトを示す）や末端のSTOPコドンが含まれていないことを確認する。
  • アライメントにイントロンやプロモーター領域など、既存のコドン置換モデルを適用できない非コード領域のヌクレオチド配列があってはならない。
  • コーディングヌクレオチド配列を直接アライメントすると、アライメントプログラムは配列のコーディング特性を考慮しないことが多いため、フレームシフト（つまり、3 の倍数ではない）ギャップが挿入されやすい。したがって、一般的には、Codon-MSAのようなコドン対応アライナーを使用するか、翻訳されたタンパク質配列をアラインメントしてから構成ヌクレオチドにマッピングし直すことを推奨する（ここでは後者の手順を採用）。
  • アラインメントに同一配列が含まれている場合、Datamonkeyは1つを除いてすべて破棄してから処理を進める。
  • HyPhyの命名規則に準拠する（スペースや特殊文字を含めない）
  • 配列が正しくアラインメントされていることを視覚的に検査することを推奨する。大きなミスアラインメントがあればを簡単に見つけることができる。 
  • 解析を実行するには、DataMonkey は少なくとも 3 つの相同なヌクレオチド配列の多重アライメントを必要とする。
  • dN および dS を推定するためのコドンベースの方法は、あらゆる配列アライメントに適用できるが、注意点がある。アライメントは、複数の分類群からサンプリングされた単一の遺伝子またはタンパク質産物、または多様な集団サンプル（例：異なる個体に感染するインフルエンザAウイルス）を表す必要がある。
    比較法は、同義置換および非同義置換の相対的な割合を推定するため、信頼性の高い推論には相当な配列多様性が必要になる。Yangらは、系統樹の全長はコドン部位あたり少なくとも1個の置換が予想される必要があると示唆しているが、望ましい配列分岐の一般的な妥当な範囲を示すことは不可能である。一方で、配列があまりにも分岐すると飽和状態、つまり枝の長さや置換パラメータを確実に推論できなくなる可能性がある。アライメントに含まれる配列数も重要で、サイトが少なすぎると、意味のある推論を行うには情報が少なすぎ、多すぎると実行に時間がかかりすぎる。他の DataMonkey ユーザーから収集された典型的なデータセットのサイズについて使用状況統計で確認できる。
  • BUSTEDでは遺伝子全体が進化しているか調べる｡aBSRELで正の選択が起きている枝が一箇所あるだけの系統樹だと､系統樹のすべての枝を平均化したとき正の選択（ω> 1）のシグナルは消えるので､そのようなケースではBUSTEDでは有意にならなかったりする。しかしBUSTEDは個別の枝ごとにテストするaBSRELが配列数が少なかったり分岐が短かったりで有意にならない少なめのデータセットでも、全枝をまとめてテストすることで統計的な検出力を確保する｡したがって､大規模なゲノムデータから、適応進化の可能性がある遺伝子を効率よくスクリーニングするなどの用途にも向いている｡

引用

Less is more: an adaptive branch-site random effects model for efficient detection of episodic diversifying selection

Martin D Smith, Joel O Wertheim, Steven Weaver, Ben Murrell, Konrad Scheffler, Sergei L Kosakovsky Pond

Mol Biol Evol. 2015 May;32(5):1342-53. 

 

Datamonkey 2.0: A Modern Web Application for Characterizing Selective and Other Evolutionary Processes
Steven Weaver, Stephen D Shank, Stephanie J Spielman, Michael Li, Spencer V Muse, Sergei L Kosakovsky Pond

Mol Biol Evol. 2018 Mar 1;35(3):773-777.

 

HyPhy 2.5-A Customizable Platform for Evolutionary Hypothesis Testing Using Phylogenies
Sergei L Kosakovsky Pond, Art F Y Poon, Ryan Velazquez, Steven Weaver, N Lance Hepler, Ben Murrell, Stephen D Shank, Brittany Rife Magalis, Dave Bouvier, Anton Nekrutenko, Sadie Wisotsky, Stephanie J Spielman, Simon D W Frost, Spencer V Muse

Mol Biol Evol. 2020 Jan 1;37(1):295-299. 

 

BUSTED

Gene-Wide Identification of Episodic Selection Free
Ben Murrell , Steven Weaver , Martin D. Smith , Joel O. Wertheim , Sasha Murrell , Anthony Aylward , Kemal Eren , Tristan Pollner , Darren P. Martin , Davey M. Smith , Konrad Scheffler , Sergei L. Kosakovsky Pond

Molecular Biology and Evolution, Volume 32, Issue 5, May 2015, Pages 1365–1371

 

参考

dN/dS 解析

 

関連