PhyloNetを使って系統推定/visualize

PhyloNetというのは系統樹の推定・作成ソフトです。このソフトはIQ-TREEなどで推定された複数の遺伝子ごとの系統樹をもとに種間のネットワークを推定することができます。が、これを使おうとして公式サイトのTutorialを見ると、インストールの次にあるのがPhyloNetの起動コマンドとinputファイルの中身が貼られており、あとはinputファイルの構造について説明があるのみで面食らいます。食らいました。

複数の系統樹を考慮して系統ネットワークを作成する、というところがミソ、ということさえ理解しておけばどうにかなると思いますが、せっかくなのでここではRAD-seqのシーケンスデータからipyrad, IQ-TREEを用いてPhyloNetのinputファイル(NEXUS形式)を作成するまでの流れを備忘録的に書いておこうと思います。
説明はipyradで処理したデータがあることを前提として進めますが、stacksなどで処理したデータでも可能と思われます。
なんとなくで構成していますので、間違いありましたらお教えいただければ...

環境

macOS 14.2.1
ipyrad  0.9.95
IQ-TREE2 2.3.4
PhyloNet 3.8.2

インストール

公式サイトにいって書いてある通りやるだけです。
これを書いているときの最新バージョンは3.8.2。
とりあえず最新のものをDLしておけばよいと思います。.Jarファイルをてきとうなディレクトリに配置すればおっけー。

Inputファイル作成の流れ

ipyradを用いてアライメントを行った場合、step7ですべてのローカスの配列が1つにまとめられたlociファイルが出力されます。これを1ローカスごとのFASTAファイルに分割します。
この分割したFASTAファイルごとにIQ-TREEで系統の推定を行い、その結果を.treefileとして出力します。この.treefile群に記述されている系統樹を一つのNEXUSファイルにまとめ、Tutorialにあるようなinputファイルを作ります。
ではやっていきましょう。

ipyrad出力のlociファイルから切り出し

ipyradで処理されたデータがあることを前提として進めます。ipyradは1~7のステップがありますが、このうちのstep7で出力される.lociファイルを使用します。.lociファイルには各ローカスの配列情報が入っているため、この配列情報をローカスごとに切り出し、FASTAファイルとして保存します。

以下のコードを任意の名前.pyで保存してターミナル上にてpython3 任意の名前.pyなどと打って実行します。すると、locus_fastaという名前のディレクトリができ、その中にlocusごとに分割されたFASTAファイルが生成されます。

ターミナル

from pathlib import Path
input_file = "input_data.loci"  # 実際のパスに合わせる
output_dir = Path("locus_fasta")
output_dir.mkdir(exist_ok=True)
with open(input_file) as f:
    locus = []
    i = 1
    for line in f:
        if line.startswith('//'):
            if locus:
                with open(output_dir / f"locus_{i}.fasta", "w") as out:
                    for entry in locus:
                        name, seq = entry
                        out.write(f">{name}\n{seq}\n")
                i += 1
                locus = []
        elif line.strip():
            parts = line.strip().split()
            if len(parts) == 2:
                locus.append((parts[0], parts[1]))

IQ-TREEで各FASTAファイルの系統推定

cdでlocus_fastaに移動し、IQ-TREEで系統の推定をします。

ターミナル

cd <パス名/locus_fasta>

for file in locus_*.fasta; do
    iqtree2 -s "$file" -m MFP -nt AUTO -pre "${file%.fasta}"
done

#IQ-TREEの細かい設定は調べてください。すまん。
#-m MFP はmodel finder。いい感じのモデルを勝手に選んでくれってやつです。

すると、.bjonj ckp.gz .iqtree .log .mldist model.gz .treefile .uniqueseq.phyなどのファイルがoutputされますが、PhyloNetのinputファイルであるNEXUS形式のファイルを作成するにあたって使用するのは.treefileのみです。この時点で.treefileには樹長のデータが入っており、これがあるとPhyloNetは文句を言って解析をしてくれないので削除する必要があります。

以下のコードを任意の名前.pyで保存してターミナル上でpython3 任意の名前.py的な感じに実行すると樹長が削除された.treefileが入ったtree_nobranchといったディレクトリが作業ディレクトリ内に作製されます。

ターミナル

import re
from pathlib import Path
input_dir = Path("<絶対パスで指定した方がよさげ/locus_fasta>")  # treefile のあるフォルダに合わせて変更
output_dir = Path("tree_nobranch")
output_dir.mkdir(exist_ok=True)
for treefile in input_dir.glob("*.treefile"):
    with open(treefile) as f:
        tree = f.read().strip()
    # 枝長情報 :0.123 のような部分を削除
    tree_no_branch = re.sub(r":-?\d+\.?\d*(e[-+]?\d+)?", "", tree)
    # 出力ファイル名
    out_file = output_dir / treefile.name
    with open(out_file, "w") as f:
        f.write(tree_no_branch + "\n")

PhyoNetのinputファイル作成

Inputファイルの作製はIQ-TREEから出力された.treefile内の系統樹をPhyloNet公式サイトにあるようなNEXUSファイルの形にまとめることとなります。先ほど樹長を削除した.treefileの中身をまとめるわけです。

以下のコードを任意の名前.pyで保存し、python3 任意の名前.pyで実行するという流れを行うと、phylonet_input.nexusというinputファイルができます。

ターミナル

cd <パス名/tree_nobranch>

from pathlib import Path
# .treefile があるフォルダ
tree_dir = Path("パス名/tree_nobranch")  
# 出力ファイル
output_file = Path("phylonet_input.nexus")
# すべての .treefile を取得
treefiles = sorted(tree_dir.glob("*.treefile"))
with open(output_file, "w") as out:
    out.write("#NEXUS\n")
    out.write("BEGIN TREES;\n")
    for i, tf in enumerate(treefiles, 1):
        with open(tf) as f:
            newick = f.read().strip()
            # 空でなければ書き出し（空行やエラー回避）
            if newick:
                out.write(f"  Tree tree_{i} = {newick}\n")
    out.write("END;\n\n")
    out.write("BEGIN PHYLONET;\n")
    out.write("  InferNetwork_ML (tree_1")
    for i in range(2, len(treefiles)+1):
        out.write(f", tree_{i}")
    out.write(") 1 -pl 6 ;\n")
    out.write("END;\n")
print(f": {len(treefiles)} 本の系統樹を含む NEXUS ファイルを生成しました: {output_file}")

細かい中身は違えど、TutorialのInferNetwork_MLにもあるような中身のファイルになるはずです。
↓こんな感じのやつ

#NEXUS

BEGIN TREES;

Tree geneTree1 = ((C,((B,D),A)),E);
Tree geneTree2 = (B,(D,(C,(A,E))));
Tree geneTree3 = (D,(B,((C,E),A)));
Tree geneTree4 = (D,((B,E),(C,A)));

END;


BEGIN PHYLONET;

InferNetwork_ML (geneTree1,geneTree2,geneTree3,geneTree4) 1 -bl 6;

END;

BEGIN PHYLONET;ブロックでコマンドを打ちます。
手法を指定（ここではInferNetwork_ML）して、推定に使用する系統樹のリストを書きます。
1 -bl;この1はreticulation数（交雑・遺伝子流動イベントの回数）です。想定されるreticulation数を書いて解析すればよいと思います。
PhyloNetはなにも指定しなかった場合コア数1で計算をするので、-bl 6など指定してコア数を指定するとよいです。

InferNetwork_ML (geneTree1,geneTree2,geneTree3,geneTree4) 1 -bl;ここの(geneTree1,geneTree2,geneTree3,geneTree4)部分にスペースを入れるとエラーを吐くので注意。

PhyloNetの実行

めでたくinputファイルができたところでPhyloNetの実行です。
特にむずかしいことはなく、PhylonNet.jarのパス inputファイルのパス outputファイルのパスといった順序で並べるだけといった雰囲気。
以下が実行のコードです。

ターミナル

java -jar パス名/PhyloNetvx_y_z.jar パス名/phylonet_input.nexus > パス名/phylonet_output.txt

みたいな感じです。
ここはTutorialを見てもわかると思います。実行すると指定したパスにoutputファイルができるはず。

visualize(図示)

無事にoutputファイルができたら、もう一仕事必要です。
outputファイルは以下のような内容と思われますが、このうちのinferred Network:のブロックにある系統樹を切り取って1行のみの.txtファイルを作ります。

InferNetwork_ML_(tree_1, tree_2, tree_3) 1 -pl 6
Results after run #1
-123.456:((((A,B),C),D),E);
Running Time (min): 1.2345
============================

Results after run #2
-123.456:((((A,B),C),D),E);
Running Time (min): 1.2345
============================

Inferred Network #1:
((A:0.2,((B:0.1)#H1:0.3):0.3,(#H1:0.2,C:0.3):0.2):0.4,D:0.5,E:0.6); 
Total log probability: -12.3456789Running time: 123.456789

この中の

((A:0.2,((B:0.1)#H1:0.3):0.3,(#H1:0.2,C:0.3):0.2):0.4,D:0.5,E:0.6); 

ここだけコピーして新たな.txtファイルを作ります。
#H1の部分が交雑(hybrid)の表記のようです。

図示にあたってはDendroscopeというソフトを使う方法と icytree ブラウザ上で動かせるツールを使う方法があります。Dendroscopeを使う場合は公式サイトに飛んでインストールしてください。

Dendroscopeはインストールが必要ですが、icytreeはオンライン上で使えるらしいです。

なお、ここで系統樹内の#Hの後に<NUMBER>::<NUMBER>といったような2重コロンが確認されるのですが、どうもこれがあるとDendroscopeが読み込んでくれません。公式Tutorialの ”Visualizing a Phylogenetic Network” を見ると、::<NUMBER>の部分をしゅどうで消去し、1つ目の樹長のみを残すといいよ、というようなことが書いてあるのでそうしましょう。

この数値は "inheritance probabilities(遺伝確率？)"らしいです。
手動で消去するか、オプション-diを使用してくれ、といったことが書かれています。
icytreeの場合は認識してくれたり認識してくれなかったりするようです。
なお、樹長の方は双方問題ないとのこと。

あとは系統樹のみ抜き出した.txtファイルをDendroscopeで開くとFigが出ます。
細かい表示の設定はまあ、てきとうに……。

2025.06.02 追記

PhyloNetの解析はざっくり触っている感じかなり遅いです。
ところで、公式のサイトを見るとPhyloNetPyなるものがあります。どうもこれはPhyloNetと比べてかなり速いようなのですが、現状はサイト内にあるTutorialに存在する解析しかできないらしい。つまり、ML法はできなさそう＆導入がかなりやっかいそうなので私自身は触れていません。
ちゃんと調べてみたら使えるよ！など、いい方法など知っている方がいましたら教えてください……。

引用文献・サイト

・PhyloNetのTutorial：https://phylogenomics.rice.edu/html/phylonetOverview.html
・Dendroscopeのマニュアル（リンク先PDF）：https://software-ab.cs.uni-tuebingen.de/download/dendroscope/manual.pdf