用於計算分子生物學的 Python 工具

Biopython 版本 1.84
© 2024。保留所有權利。
在 GitHub 上編輯此頁面

Bio.Phylo 食譜。

這裡有一些使用 Bio.Phylo 來完成一些常見任務的範例。這些函數中的一些可能會在之後的版本中加入 Biopython，但是您可以在您的程式碼中使用 Biopython 1.54。

方便函數

取得演化支的父節點

Bio.Phylo 中的 Tree 資料結構不會儲存每個演化支的父節點參考。相反地，可以使用 get_path 方法來追蹤從樹根到所選演化支的父子連結路徑。

def get_parent(tree, child_clade):
    node_path = tree.get_path(child_clade)
    return node_path[-2]


# Select a clade
myclade = tree.find_clades("foo").next()
# Test the function
parent = get_parent(tree, myclade)
assert myclade in parent

請注意，get_path 的執行時間與樹的大小呈線性關係 - 也就是說，為了獲得最佳效能，請勿在時間敏感的迴圈內呼叫 get_parent 或 get_path。如果可以，請在迴圈外部呼叫 get_path，並在該函數傳回的列表中查找父節點。

或者，如果您需要重複查找任意樹元素的父節點，請建立一個字典，將所有節點映射到它們的父節點。

def all_parents(tree):
    parents = {}
    for clade in tree.find_clades(order="level"):
        for child in clade:
            parents[child] = clade
    return parents


# Example
parents = all_parents(tree)
myclade = tree.find_clades("foo").next()
parent_of_myclade = parents[myclade]
assert myclade in parent_of_myclade

依名稱索引演化支

對於大型樹狀結構，能夠依名稱或其他唯一識別碼選擇演化支可能會很有用，而不是在每次操作期間搜尋整棵樹。

def lookup_by_names(tree):
    names = {}
    for clade in tree.find_clades():
        if clade.name:
            if clade.name in names:
                raise ValueError("Duplicate key: %s" % clade.name)
            names[clade.name] = clade
    return names

現在您可以透過名稱在常數時間內檢索演化支

tree = Phylo.read("ncbi_taxonomy.xml", "phyloxml")
names = lookup_by_names(tree)
for phylum in ("Apicomplexa", "Euglenozoa", "Fungi"):
    print("Phylum size: %d" % len(names[phylum].get_terminals()))

一個潛在的問題：上述 lookup_by_names 的實作不包括未命名的演化支，通常是內部節點。我們可以透過為每個演化支添加一個唯一識別碼來解決這個問題。在這裡，所有演化支名稱都加上一個唯一的數字前綴（這對於搜尋也很有用）

def tabulate_names(tree):
    names = {}
    for idx, clade in enumerate(tree.find_clades()):
        if clade.name:
            clade.name = "%d_%s" % (idx, clade.name)
        else:
            clade.name = str(idx)
        names[clade.name] = clade
    return names

計算相鄰終端之間的距離

由 Joel Berendzen 建議

import itertools


def terminal_neighbor_dists(self):
    """Return a list of distances between adjacent terminals."""

    def generate_pairs(self):
        pairs = itertools.tee(self)
        next(pairs[1])
        return zip(pairs[0], pairs[1])

    return [self.distance(*i) for i in generate_pairs(self.find_clades(terminal=True))]

測試「半終端」演化支

由 Joel Berendzen 建議

現有的樹方法 is_preterminal 在所有直接後代都是終端時傳回 True。這個程式碼片段會在任何直接後代是終端時傳回 True，但如果給定的演化支本身是終端時仍然會傳回 False。

def is_semipreterminal(clade):
    """True if any direct descendent is terminal."""
    for child in clade:
        if child.is_terminal():
            return True
    return False

在 Python 2.5 和更新的版本中，可以使用內建的 any 函數來簡化

def is_semipreterminal(clade):
    return any(child.is_terminal() for child in clade)

比較樹狀結構

待辦事項

共識方法

待辦事項

定根方法

Tree 類別（而非 TreeMixin）上的基本方法是 root_with_outgroup

tree = Phylo.read("example.nwk", "newick")
print(tree)
# ...
tree.root_with_outgroup({"name": "A"})  # Operates in-place
print(tree)

通常您會希望外群是一個單系群，而不是單一分類群。這不會自動檢查，但您可以使用 is_monophyletic 方法自行檢查。

為了節省一些打字時間，請嘗試將查詢保留在單獨的列表中並重複使用

outgroup = [{"name": taxon_name} for taxon_name in ("E", "F", "G")]
if tree.is_monophyletic(outgroup):
    tree.root_with_outgroup(*outgroup)
else:
    raise ValueError("outgroup is paraphyletic")

待辦事項

圖形

待辦事項

匯出到其他類型

透過 Rpy2 轉換為 'ape' 樹

R 統計程式設計環境透過 ape 套件以及其他幾個建立在 ape 之上的套件，提供對系統發育的支持。Python 套件 rpy2 提供 R 和 Python 之間的介面，因此可以將 Bio.Phylo 樹轉換為 ape 樹物件。

import tempfile
from rpy2.robjects import r


def to_ape(tree):
    """Convert a tree to the type used by the R package `ape`, via rpy2.

    Requirements:
        - Python package `rpy2`
        - R package `ape`
    """
    with tempfile.NamedTemporaryFile() as tmpf:
        Phylo.write(tree, tmpf, "newick")
        tmpf.flush()
        rtree = r(
            """
            library('ape')
            read.tree('%s')
            """
            % tmpf.name
        )
    return rtree

確認它有效

>>> from StringIO import StringIO
>>> from Bio import Phylo
>>> tree = Phylo.read(StringIO("(A,(B,C),(D,E));"), "newick")
>>> rtree = to_ape(tree)
>>> len(rtree)
3
>>> print(r.summary(rtree))
Phylogenetic tree: structure(list(edge = structure(c(6, 6, 7, 7, 6, 8, 8, 1, 7,  2, 3, 8, 4, 5),
 .Dim = c(7L, 2L)), tip.label = c("A", "B", "C",  "D", "E"), Nnode = 3L),
 .Names = c("edge", "tip.label", "Nnode" ), class="phylo")

  Number of tips: 5
  Number of nodes: 3
  No branch lengths.
  No root edge.
  Tip labels: A
              B
              C
              D
              E
  No node labels.
NULL
>>> r.plot(rtree)

請參閱 rpy2 文件以獲得進一步的指導。

轉換為 DendroPy 或 PyCogent 樹

Biopython、DendroPy 和 PyCogent 使用的樹物件不同。儘管如此，所有三個工具包都支援 Newick 檔案格式，因此，透過使用一個庫寫入暫存檔或 StringIO 物件，然後使用另一個庫再次讀取相同的字串，可以在該層級實現互通性。

from Bio import Phylo
import cogent

Phylo.write(bptree, "mytree.nwk", "newick")  # Biopython tree
ctree = cogent.LoadTree("mytree.nwk")  # PyCogent tree

import dendropy

# Create or load a tree in DendroPy
dtree = dendropy.Tree.get_from_string("(A, (B, C), (D, E))", "newick")
dtree.write_to_path("tmp.nwk", "newick", suppress_rooting=True)
# Load the same tree in Biopython
bptree = Phylo.read("tmp.nwk", "newick")

轉換為 NumPy 陣列或矩陣

鄰接矩陣： 如果存在父子關係，則單元格為 1 (true)，否則為 0 (false)。

import numpy as np


def to_adjacency_matrix(tree):
    """Create an adjacency matrix (NumPy array) from clades/branches in tree.

    Also returns a list of all clades in tree ("allclades"), where the position
    of each clade in the list corresponds to a row and column of the numpy
    array: a cell (i,j) in the array is 1 if there is a branch from allclades[i]
    to allclades[j], otherwise 0.

    Returns a tuple of (allclades, adjacency_matrix) where allclades is a list
    of clades and adjacency_matrix is a NumPy 2D array.
    """
    allclades = list(tree.find_clades(order="level"))
    lookup = {}
    for i, elem in enumerate(allclades):
        lookup[elem] = i
    adjmat = np.zeros((len(allclades), len(allclades)))
    for parent in tree.find_clades(terminal=False, order="level"):
        for child in parent.clades:
            adjmat[lookup[parent], lookup[child]] = 1
    if not tree.rooted:
        # Branches can go from "child" to "parent" in unrooted trees
        adjmat = adjmat + adjmat.transpose()
    return (allclades, np.matrix(adjmat))

距離矩陣： 如果存在分支，則單元格值為分支長度，否則為無窮大（這與圖形演算法配合良好）。

import numpy as np


def to_distance_matrix(tree):
    """Create a distance matrix (NumPy array) from clades/branches in tree.

    A cell (i,j) in the array is the length of the branch between allclades[i]
    and allclades[j], if a branch exists, otherwise infinity.

    Returns a tuple of (allclades, distance_matrix) where allclades is a list of
    clades and distance_matrix is a NumPy 2D array.
    """
    allclades = list(tree.find_clades(order="level"))
    lookup = {}
    for i, elem in enumerate(allclades):
        lookup[elem] = i
    distmat = np.repeat(np.inf, len(allclades) ** 2)
    distmat.shape = (len(allclades), len(allclades))
    for parent in tree.find_clades(terminal=False, order="level"):
        for child in parent.clades:
            if child.branch_length:
                distmat[lookup[parent], lookup[child]] = child.branch_length
    if not tree.rooted:
        distmat += distmat.transpose
    return (allclades, np.matrix(distmat))

增強功能

待辦事項