チュートリアル

SplitAligner を、はじめから終わりまで

種の系統樹と遺伝子系統樹の集合を枝座標行列に変換する。各枝は数値、または明示的な欠測状態——構造的・融合関連・トポロジー由来・残余 NA——を持つ。本チュートリアルは同梱のトイ例を使うので、すべての出力を自分の環境で再現できる。

SplitAligner は小さく依存の少ない Perl プログラムである。コンパイルは不要——リポジトリをクローンして SplitAligner.pl を実行するだけ。以下のコマンドは examples/run.sh の処理手順に対応している。
考え方

枝とは、それが作るスプリットである

SplitAligner は枝を、描かれた位置や節点の順序ではなく、それが分類群に誘導する二分割(スプリット)で同定する。一部の分類群 Tg しか含まない遺伝子 g について、種の系統樹の各枝 b を、実際に存在する分類群へ射影する:

σg(b) = ( Ab ∩ Tg | Bb ∩ Tg )

種樹の制限木とトポロジーが一致する基準層では、削除されていない二本の原始枝の射影スプリットが一致することと、同じ簡約辺へ写ることは同値である。これは元の二枝が同じ枝になったという意味ではない。その遺伝子の分類群集合では二枝を独立に区別できなくなった、つまり融合したという意味である。射影後の片側が空なら、代わりに構造的欠測となる。経験的な free トポロジー遺伝子樹については、独立に観測可能な各射影スプリットが回収されたかを別に判定する。

表示上の根は座標ではない。SplitAligner の生物学的座標は正準無根スプリットであり、描画された根の位置、内部節点の番号、走査順序、子の順序ではない。二つの Newick 表現が同一の重み付き無根スプリット対象を符号化しているならば、表示上の根や直列化だけを変えても、基礎にある生物学的なスプリットの同一性は変わらない。

可読性のための B エイリアスは直列化に固有であり、独立に直列化された種樹の間では変わりうる。独立した表現どうしを比較する場合、権威ある同一性は B 番号そのものではなく正準スプリット鍵である。

根の位置は、方向や祖先関係を要する解析では依然として生物学的な意味をもちうる。それが SplitAligner の枝座標の同一性を定義しないというだけである。

表現 A — 内部辺で根を置く AB CD 表現 B — 枝 A で根を置く AB CD 正準スプリット鍵は同一 A,B | C,D B エイリアスは異なりうる

ワークフローには二つのモードがある。matrix が座標系を構築し、全遺伝子をそこへ射影する。finalize が、対応する証拠に基づいて欠測状態を分類する。

1

インストール

Perl 5 が必要である(macOS や多くの Linux には最初から入っている)。SplitAligner はコアモジュールのみを使い、CPAN パッケージのインストールは不要である。

# リポジトリをクローン
git clone https://github.com/wujiaqi06/SplitAligner
cd SplitAligner

# 動作確認(使い方バナーが表示される)
perl SplitAligner.pl

インストールはこれで全部である。メインの制御は SplitAligner.pl。各ステップの補助スクリプトは scripts/ にあり、自動的に呼び出される。

2

入力を用意する

種の系統樹 — Newick 一本。枝長やサポート値は許容され、スプリットに基づくマッピングでは単純に無視される。

((A:0.1,B:0.2):0.2,(C:0.1,D:0.1):0.1):0.1;

遺伝子系統樹 — Newick、一行一レコード。各行は遺伝子識別子で始まり、その直後に系統樹が続く。分類群名は種の系統樹と一致していなければならない。

GeneA((A:0.1,B:0.2):0.2,(C:0.1,D:0.1):0.1):0.1;
GeneB((A:0.2,C:0.1):0.1,(B:0.1,D:0.2):0.1):0.1;

SplitAligner は二種類の遺伝子系統樹入力を区別し、通常は両方を与える:

  • free — 各遺伝子を独立に推定した系統樹(トポロジーが種樹と食い違ってよい);
  • fix — 種樹トポロジーに拘束して推定した系統樹(枝長のみ)。

この二つを比較することが、のちにカバレッジ由来の欠測と不一致由来の欠測を分ける鍵になる。

3

クイックスタート:トイ例

リポジトリには examples/302mammal/ に小さな例(free / fix 遺伝子系統樹のサブセット)が同梱されており、端から端までの動作確認を素早く行える。リポジトリのルートから:

bash examples/run.sh toy

これは以下の三回のプログラム呼び出しを順に実行する。本チュートリアルの残りでは、各ステップが何をし何を生成するかを解きほぐす。

4

matrix モード — 座標系を構築する

matrix は遺伝子系統樹の集合ごとに一度実行する。種樹の枝にラベルを付け、全遺伝子をそこへ射影し、遺伝子 × 枝の行列を書き出す。

# リポジトリのルートから開始
cd examples/302mammal

perl ../../SplitAligner.pl --mode matrix \
  --species input/speciesTree302.nwk \
  --gene input/free_tree.examples.nwk \
  --label free

perl ../../SplitAligner.pl --mode matrix \
  --species input/speciesTree302.nwk \
  --gene input/fix_tree.examples.nwk \
  --label fix

主な出力(--label が接頭辞):

  • species_tree.splits.txt — 標準的な枝座標系(種樹バックボーン上の B1, B2, …)。全遺伝子が測られる共通の軸である。
  • <label>.matrix_no_fuse.txt — 原始枝のみ。
  • <label>.matrix_with_fuse.txt — それに加えて、分類群の剪定で融合する枝の明示的な融合枝カラム(例 B12|B47)。
  • <label>_splits/<label>_split_branch_label/ — 遺伝子ごとのスプリット表現と、そのバックボーンへの対応付け。

各セルは、枝の射影スプリットが回収された場合はその枝長(数値)を、そうでない場合は NA を——まだ——持つ。この NA に名前を与えるのが次のステップである。

gene   B1           B2   B3           B4           ...
A1BG   0.1809315089 NA   0.0047384154 0.0057110032 ...
A1CF   0.0376651404 0.0583312140 0.0032591533 ...
free.matrix_with_fuse.txt(一部抜粋)
5

finalize モード — 欠測状態を分類する

finalize に両方の行列を渡す。まず、数値を持つ融合座標で説明できる原始枝の NANA_fuse として識別し、次に共有遺伝子上で free と fix の行列を比較して、構造的欠測(NA_struct)とトポロジー由来の欠測(NA_topo)を切り分ける。fix 側に判定可能な数値的証拠がなければ、残余 NA を保持する。--species_tree を渡すと枝ごとのサポートも計算する。

perl ../../SplitAligner.pl --mode finalize \
  --free free.matrix_with_fuse.txt \
  --fix  fix.matrix_with_fuse.txt \
  --final_label final \
  --species_tree species_tree.forSplit.nwk

出力:

  • final.free.na_classified.txt / final.fix.na_classified.txt — 状態分類済みの行列。fix 側に数値的証拠がなければ、残余 NA が保持される;
  • <input>.na_fuse.txt — 融合救済の中間テーブル;
  • final.support_b.txtspecies_tree.support_b.nwk — 枝ごとのサポート(ステップ 7)。

比較は両方の入力に存在する遺伝子についてのみ定義される。共有遺伝子が無ければ、SplitAligner は推測せずエラーで停止する。

gene   B1           B2         B3           ...
A1BG   0.1809315089 NA_struct  0.0047384154 ...
final.free.na_classified.txt — B2 の一般的な NA に名前が付いた
6

状態の読み方

すべての遺伝子 × 枝のセルは、最終的に以下のいずれか一つの値を取る。残余 NA は意図的な記帳状態である。fix 側の数値的証拠なしに NA_topo へ格上げすることはない。

セルの値意味
数値Mapped(マップ済)。枝の射影スプリットがこの遺伝子で観測された。値はその枝長。
NA_fuse分類群の剪定後、この原始枝は独立には観測できない。融合座標は数値の枝長を持ち得るが、その値は各原始枝の独立した推定値ではない。
NA_struct射影された側が空になり、この遺伝子について枝が射影上の同一性を持たず評価不能(カバレッジの効果)。
NA_topofix 側には原始枝の数値的証拠があるが、その射影スプリットが推定された free トポロジー樹では回収されない——トポロジー由来の不一致と整合する状態。
NA残余:fix 側に判定できる数値的根拠が無いときに保持される。

これらのカテゴリは軸を分割します——台帳全体を足し合わせる:

|C| = Mapped + NA_fuse + NA_struct + NA_topo + 残余 NA

NA理由ごとに分ければ、下流の速度・選択解析で各状態を明示的に処理、フィルタリング、監査できる。すべての欠測を一つの不透明な値として扱う必要はない。

7

枝ごとの Support

--species_tree を渡すと、finalizefinal.support_b.txt を書き出す。バックボーンの各枝について、共有遺伝子上で fix の系統樹に対する free の系統樹での数値的根拠の保持度を示す。

branch_id  branch_type  n_shared_genes  n_fix_non_na  n_free_non_na  support_percent  discordance_percent
B1         terminal     5               5             5              100.0000000000   0.0000000000
B2         terminal     5               4             4              100.0000000000   0.0000000000
B3         terminal     5               4             4              100.0000000000   0.0000000000
final.support_b.txt(先頭)

対になる species_tree.support_b.nwk は、各内部枝の Support をブートストラップ位置に書き込んだ標準 Newick 樹である。任意の樹ビューア(例 FigTree)で開けば、一致度をトポロジー上に直接投影して見られる。

8

スケールアップ

同じ二つのモードが本番規模でも動く。同梱のプレプリント例は、302 種の哺乳類にわたる 2,275 本の free トポロジーと 2,275 本の fix トポロジーの遺伝子系統樹を使う:

bash examples/run.sh preprint

共通の一致性サブ問題において、SplitAligner の Support(b)PhyParts と高い一致を示す(299 本の内部枝で Pearson r = 0.99997)。論文で報告した同一マシン上の比較では、SplitAligner の全工程は 9 分 52 秒、PhyParts の quick concordance 設定は 1 時間 50 分 53 秒で、このデータセットでは約 11 倍高速であった。これは校正と実行時間の比較であり、機能的な同等性を主張するものではない。

9

Catnip10 を動かす

遺伝子樹の不一致を加えず、分類群の剪定だけを作用させる。完全な B1–B17 軸は固定されたまま、削除のたびにグラフで定まる台帳が observedNA_fuseNA_struct へ更新される。この実験では NA_topo は常にゼロである。