# ----------------------------------------------------------------------------
# Copyright (c) 2013--, scikit-bio development team.
#
# Distributed under the terms of the Modified BSD License.
#
# The full license is in the file LICENSE.txt, distributed with this software.
# ----------------------------------------------------------------------------

from warnings import warn, simplefilter
from operator import or_, itemgetter
from copy import copy, deepcopy
from itertools import combinations
from functools import reduce
from collections import defaultdict

import numpy as np
import pandas as pd
from scipy.spatial.distance import correlation

from skbio._base import SkbioObject
from skbio.stats.distance import DistanceMatrix
from skbio.tree._exception import (
    NoLengthError,
    DuplicateNodeError,
    NoParentError,
    MissingNodeError,
    TreeError,
)
from skbio.util import RepresentationWarning
from skbio.util._decorator import classonlymethod
from skbio.util._warning import _warn_deprecated


def distance_from_r(m1, m2):
    r"""Estimate distance as (1-r)/2: neg correl = max distance.

    Parameters
    ----------
    m1 : DistanceMatrix
        a distance matrix to compare
    m2 : DistanceMatrix
        a distance matrix to compare

    Returns
    -------
    float
        The distance between m1 and m2

    """
    return correlation(m1.data.flat, m2.data.flat) / 2


class TreeNode(SkbioObject):
    r"""Representation of a node within a tree.

    A `TreeNode` instance stores links to its parent and optional children
    nodes. In addition, the `TreeNode` can represent a `length` (e.g., a
    branch length) between itself and its parent. Within this object, the use
    of "children" and "descendants" is frequent in the documentation. A child
    is a direct descendant of a node, while descendants are all nodes that are
    below a given node (e.g., grand-children, etc).

    Parameters
    ----------
    name : str or None
        A node can have a name. It is common for tips in particular to have
        names, for instance, in a phylogenetic tree where the tips correspond
        to species.
    length : float, int, or None
        Length of the branch connecting this node to its parent. Can represent
        ellapsed time, amount of mutations, or other measures of evolutionary
        distance.
    support : float, int, or None
        Support value of the branch connecting this node to its parent. Can be
        bootstrap value, posterior probability, or other metrics measuring the
        confidence or frequency of this branch.
    parent : TreeNode or None
        Connect this node to a parent
    children : list of TreeNode or None
        Connect this node to existing children

    """

    default_write_format = "newick"
    _exclude_from_copy = set(["parent", "children", "_tip_cache", "_non_tip_cache"])

    def __init__(
        self, name=None, length=None, support=None, parent=None, children=None
    ):
        self.name = name
        self.length = length
        self.support = support
        self.parent = parent
        self.children = []
        self.id = None
        if children is not None:
            self.extend(children)

    def __repr__(self):
        r"""Return summary of the tree.

        Returns
        -------
        str
            A summary of this node and all descendants

        Notes
        -----
        This method returns the name of the node and a count of tips and the
        number of internal nodes in the tree

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c, d)root;"])
        >>> repr(tree)
        '<TreeNode, name: root, internal node count: 1, tips count: 3>'

        """
        nodes = [n for n in self.traverse(include_self=False)]
        n_tips = sum([n.is_tip() for n in nodes])
        n_nontips = len(nodes) - n_tips
        classname = self.__class__.__name__
        name = self.name if self.name is not None else "unnamed"

        return "<%s, name: %s, internal node count: %d, tips count: %d>" % (
            classname,
            name,
            n_nontips,
            n_tips,
        )

    def __str__(self):
        r"""Return string version of self, with names and distances.

        Returns
        -------
        str
            Returns a Newick representation of the tree

        See Also
        --------
        read
        write

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c);"])
        >>> str(tree)
        '((a,b)c);\n'

        """
        return str("".join(self.write([])))

    def __iter__(self):
        r"""Node iter iterates over the `children`."""
        return iter(self.children)

    def __len__(self):
        return len(self.children)

    def __getitem__(self, i):
        r"""Node delegates slicing to `children`."""
        return self.children[i]

    def _adopt(self, node):
        r"""Update `parent` references but does NOT update `children`."""
        if node.parent is not None:
            node.parent.remove(node)
        node.parent = self
        return node

    def append(self, node):
        r"""Add a node to self's children.

        Parameters
        ----------
        node : TreeNode
            Node to add as a child.

        See Also
        --------
        extend

        Notes
        -----
        ``append`` will invalidate any node lookup caches, remove the node's
        parent if it exists, set the parent of node to self, and add the node
        to self's children.

        Examples
        --------
        >>> from skbio import TreeNode
        >>> root = TreeNode(name="root")
        >>> child1 = TreeNode(name="child1")
        >>> child2 = TreeNode(name="child2")
        >>> root.append(child1)
        >>> root.append(child2)
        >>> print(root)
        (child1,child2)root;
        <BLANKLINE>

        """
        self.invalidate_caches()
        self.children.append(self._adopt(node))

    def extend(self, nodes):
        r"""Add a list of nodes to self's children.

        Parameters
        ----------
        nodes : list of TreeNode
            Nodes to add as children.

        See Also
        --------
        append

        Notes
        -----
        ``extend`` will invalidate any node lookup caches, remove existing
        parents of the nodes if they have any, set their parents to self
        and add the nodes to the children of self.

        Examples
        --------
        >>> from skbio import TreeNode
        >>> root = TreeNode(name="root")
        >>> root.extend([TreeNode(name="child1"), TreeNode(name="child2")])
        >>> print(root)
        (child1,child2)root;
        <BLANKLINE>

        """
        self.invalidate_caches()
        self.children.extend([self._adopt(n) for n in nodes[:]])

    def insert(self, node, distance=None, branch_attrs=[]):
        r"""Insert a node into the branch connecting self and its parent.

        .. versionadded:: 0.6.2

        Parameters
        ----------
        node : TreeNode
            Node to insert.
        distance : float, int or None, optional
            Distance between self and the insertion point. Must not exceed
            ``self.length``. If ``None`` whereas ``self.length`` is not
            ``None``, will insert at the midpoint of the branch.
        branch_attrs : iterable of str, optional
            Attributes of self that should be transferred to the inserted node
            as they are considered as attributes of the branch. ``support``
            will be automatically included as it is always a branch attribute.

        Raises
        ------
        NoParentError
            If self has no parent.
        ValueError
            If distance is specified but branch has no length.
        ValueError
            If distance exceeds branch length.

        See Also
        --------
        append

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a:1,b:2)c:4,d:5)e;"])
        >>> print(tree.ascii_art())
                            /-a
                  /c-------|
        -e-------|          \-b
                 |
                  \-d

        >>> tree.find("c").insert(TreeNode("x"))
        >>> print(tree.ascii_art())
                                      /-a
                  /x------- /c-------|
        -e-------|                    \-b
                 |
                  \-d
        >>> tree.find("c").length
        2.0
        >>> tree.find("x").length
        2.0

        """
        if (parent := self.parent) is None:
            raise NoParentError("Self has no parent.")

        # detach node from original tree if applicable
        if node.parent is not None:
            node.parent.remove(node)

        # See also `_adopt`. The current code replaces the node at the same
        # position in the parent's list of children, instead of appending to
        # the end. Additionally, the current code performs cache invalidation
        # only once.
        self.invalidate_caches()

        # replace self with node in the parent's list of children
        node.parent = parent
        for i, curr_node in enumerate(parent.children):
            if curr_node is self:
                parent.children[i] = node

        # add self to the beginning of the node's list of children
        self.parent = node
        node.children.insert(0, self)

        # transfer branch attributes to new node
        branch_attrs = set(branch_attrs)
        branch_attrs.add("support")
        branch_attrs.discard("length")
        for attr in branch_attrs:
            setattr(node, attr, getattr(self, attr, None))

        # determine insertion point
        if distance is None:
            if self.length is None:
                node.length = None
            else:
                self.length *= 0.5
                node.length = self.length
        else:
            if self.length is None:
                raise ValueError("Distance is provided but branch has no length.")
            elif distance > self.length:
                raise ValueError("Distance cannot exceed branch length.")
            node.length = self.length - distance
            self.length = distance

    def pop(self, index=-1):
        r"""Remove and return a child node by its index position from self.

        Parameters
        ----------
        index : int
            The index position in ``children`` to pop.

        Returns
        -------
        TreeNode
            The popped child node.

        See Also
        --------
        remove
        remove_deleted

        Notes
        -----
        All node lookup caches are invalidated, and the parent reference for
        the popped node will be set to ``None``.

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["(a,b)c;"])
        >>> print(tree.pop(0))
        a;
        <BLANKLINE>

        """
        return self._remove_node(index)

    def _remove_node(self, idx):
        r"""Perform node removal.

        The actual (and only) method that performs node removal.

        """
        self.invalidate_caches()
        node = self.children.pop(idx)
        node.parent = None
        return node

    def remove(self, node):
        r"""Remove a node from self.

        Remove a `node` from `self` by identity of the node.

        Parameters
        ----------
        node : TreeNode
            The node to remove from self's children

        Returns
        -------
        bool
            `True` if the node was removed, `False` otherwise

        See Also
        --------
        pop
        remove_deleted

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["(a,b)c;"])
        >>> tree.remove(tree.children[0])
        True

        """
        for i, curr_node in enumerate(self.children):
            if curr_node is node:
                self._remove_node(i)
                return True
        return False

    def remove_deleted(self, func):
        r"""Delete nodes in which `func(node)` evaluates `True`.

        Remove all descendants from `self` that evaluate `True` from `func`.
        This has the potential to drop clades.

        Parameters
        ----------
        func : a function
            A function that evaluates `True` when a node should be deleted

        See Also
        --------
        pop
        remove

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["(a,b)c;"])
        >>> tree.remove_deleted(lambda x: x.name == 'b')
        >>> print(tree)
        (a)c;
        <BLANKLINE>

        """
        for node in self.traverse(include_self=False):
            if func(node):
                node.parent.remove(node)

    def prune(self):
        r"""Reconstruct correct topology after nodes have been removed.

        Internal nodes with only one child will be removed and new connections
        will be made to reflect change. This method is useful to call
        following node removals as it will clean up nodes with singular
        children.

        Names and properties of singular children will override the names and
        properties of their parents following the prune.

        Node lookup caches are invalidated.

        See Also
        --------
        shear
        remove
        pop
        remove_deleted

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c,(d,e)f)root;"])
        >>> to_delete = tree.find('b')
        >>> tree.remove_deleted(lambda x: x == to_delete)
        >>> print(tree)
        ((a)c,(d,e)f)root;
        <BLANKLINE>
        >>> tree.prune()
        >>> print(tree)
        ((d,e)f,a)root;
        <BLANKLINE>

        """
        # build up the list of nodes to remove so the topology is not altered
        # while traversing
        nodes_to_remove = []
        for node in self.traverse(include_self=False):
            if len(node.children) == 1:
                nodes_to_remove.append(node)

        # clean up the single children nodes
        for node in nodes_to_remove:
            child = node.children[0]

            if child.length is None or node.length is None:
                child.length = child.length or node.length
            else:
                child.length += node.length

            if node.parent is None:
                continue

            node.parent.append(child)
            node.parent.remove(node)

        # if a single descendent from the root, the root adopts the childs
        # properties. we can't "delete" the root as that would be deleting
        # self.
        if len(self.children) == 1:
            node_to_copy = self.children[0]
            efc = self._exclude_from_copy
            for key in node_to_copy.__dict__:
                if key not in efc:
                    self.__dict__[key] = deepcopy(node_to_copy.__dict__[key])
            self.remove(node_to_copy)
            self.extend(node_to_copy.children)

    def shear(self, names):
        """Remove tips until the tree just has the desired tip names.

        Parameters
        ----------
        names : Iterable of str
            The tip names on the tree to keep

        Returns
        -------
        TreeNode
            The resulting tree

        Raises
        ------
        ValueError
            If the names do not exist in the tree

        See Also
        --------
        prune
        remove
        pop
        remove_deleted

        Examples
        --------
        >>> from skbio import TreeNode
        >>> t = TreeNode.read(['((H:1,G:1):2,(R:0.5,M:0.7):3);'])
        >>> sheared = t.shear(['G', 'M'])
        >>> print(sheared)
        (G:3.0,M:3.7);
        <BLANKLINE>

        """
        tcopy = self.copy()
        all_tips = {n.name for n in tcopy.tips()}
        ids = set(names)

        if not ids.issubset(all_tips):
            raise ValueError("ids are not a subset of the tree.")

        marked = set()
        for tip in tcopy.tips():
            if tip.name in ids:
                marked.add(tip)
                for anc in tip.ancestors():
                    if anc in marked:
                        break
                    else:
                        marked.add(anc)

        for node in list(tcopy.traverse()):
            if node not in marked:
                node.parent.remove(node)

        tcopy.prune()

        return tcopy

    def _copy(self, deep, memo):
        """Return a copy of self."""
        _copy = deepcopy if deep else copy
        _args = [memo] if deep else []

        def __copy_node(node_to_copy):
            """Copy a node."""
            # this is _possibly_ dangerous, we're assuming the node to copy is
            # of the same class as self, and has the same exclusion criteria.
            # however, it is potentially dangerous to mix TreeNode subclasses
            # within a tree, so...
            result = self.__class__()
            efc = self._exclude_from_copy
            for key in node_to_copy.__dict__:
                if key not in efc:
                    result.__dict__[key] = _copy(node_to_copy.__dict__[key], *_args)
            return result

        root = __copy_node(self)
        nodes_stack = [[root, self, len(self.children)]]

        while nodes_stack:
            # check the top node, any children left unvisited?
            top = nodes_stack[-1]
            new_top_node, old_top_node, unvisited_children = top

            if unvisited_children:
                top[2] -= 1
                old_child = old_top_node.children[-unvisited_children]
                new_child = __copy_node(old_child)
                new_top_node.append(new_child)
                nodes_stack.append([new_child, old_child, len(old_child.children)])
            else:  # no unvisited children
                nodes_stack.pop()
        return root

    def __copy__(self):
        """Return a shallow copy."""
        return self._copy(False, {})

    def __deepcopy__(self, memo):
        """Return a deep copy."""
        return self._copy(True, memo)

    def copy(self, deep=True):
        r"""Return a copy of self using an iterative approach.

        Parameters
        ----------
        deep : bool, optional
            Whether perform a deep (``True``, default) or shallow (``False``)
            copy of node attributes.

            .. versionadded:: 0.6.2

            .. note:: The default value will be changed to ``False`` in 0.7.0.

        Returns
        -------
        TreeNode
            A new copy of self.

        See Also
        --------
        unrooted_copy

        Notes
        -----
        This method iteratively copies the current node and its descendants.
        That is, if the current node is not the root of the tree, only the
        subtree below the node, instead of the entire tree, will be copied.

        All nodes and their attributes will be copied. The copies are new
        objects rather than references to the original objects. The distinction
        between deep and shallow copies only applies to each node attribute.

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c,(d,e)f)root;"])
        >>> tree_copy = tree.copy()
        >>> tree_nodes = set([id(n) for n in tree.traverse()])
        >>> tree_copy_nodes = set([id(n) for n in tree_copy.traverse()])
        >>> print(len(tree_nodes.intersection(tree_copy_nodes)))
        0

        """
        return self._copy(deep, {})

    def deepcopy(self):
        r"""Return a deep copy of self using an iterative approach.

        Returns
        -------
        TreeNode
            A new deep copy of self.

        See Also
        --------
        copy

        Notes
        -----
        ``deepcopy`` is equivalent to ``copy`` with ``deep=True``, which is
        currently the default behavior of the latter.

        Warnings
        --------
        ``deepcopy`` is deprecated as of ``0.6.2``. Use ``copy`` instead.

        """
        msg = "Use copy instead."
        _warn_deprecated(self.__class__.deepcopy, "0.6.2", msg)

        return self._copy(True, {})

    def unrooted_copy(
        self,
        parent=None,
        branch_attrs={"name", "length", "support"},
        root_name="root",
        deep=False,
    ):
        r"""Walk the tree unrooted-style and return a copy.

        Parameters
        ----------
        parent : TreeNode or None
            Direction of walking (from parent to self). If specified, walking
            to the parent will be prohibited.

        branch_attrs : set of str, optional
            Attributes of ``TreeNode`` objects that should be considered as
            branch attributes during the operation.

            .. versionadded:: 0.6.2

            .. note:: ``name`` will be removed from the default in 0.7.0, as
                it is usually considered as an attribute of the node instead of
                the branch.

        root_name : str or None, optional
            Name for the new root node, if it doesn't have one.

            .. versionadded:: 0.6.2

            .. note:: This parameter will be removed in 0.7.0, and the root
                node will not be renamed.

        deep : bool, optional
            Whether perform a shallow (``False``, default) or deep (``True``)
            copy of node attributes.

            .. versionadded:: 0.6.2

        Returns
        -------
        TreeNode
            A new copy of the tree rooted at the given node.

            .. versionchanged:: 0.6.2

                Node attributes other than name and length will also be copied.

        Warnings
        --------
        The default behavior of ``unrooted_copy`` is subject to change in
        0.7.0. The new default behavior can be achieved by specifying
        ``branch_attrs={"length", "support"}, root_name=None``.

        See Also
        --------
        copy
        unrooted_move

        Notes
        -----
        This method recursively walks a tree from a given node in an unrooted
        style (i.e., directions of branches are not assumed), and copies each
        node it visits, such that the copy of the given node becomes the root
        node of a new tree and the copies of all other nodes are re-positioned
        accordingly, whereas the topology of the new tree will be identical to
        the existing one.

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,(b,c)d)e,(f,g)h)i;"])
        >>> new_tree = tree.find('d').unrooted_copy()
        >>> print(new_tree)
        (b,c,(a,((f,g)h)e)d)root;
        <BLANKLINE>

        """
        # future warning
        if branch_attrs == {"name", "length", "support"} and root_name == "root":
            func = self.__class__.unrooted_copy
            if not hasattr(func, "warned"):
                simplefilter("once", FutureWarning)
                warn(
                    "The default behavior of `unrooted_copy` is subject to change in "
                    "0.7.0. The new default behavior can be achieved by specifying "
                    '`branch_attrs={"length", "support"}, root_name=None`.',
                    FutureWarning,
                )
                func.warned = True

        _copy = deepcopy if deep else copy

        # identify neighbors (adjacent nodes) of self, excluding the incoming node
        neighbors = self.neighbors(ignore=parent)

        # recursively copy each neighbor; they will become outgoing nodes (children)
        children = [
            c.unrooted_copy(
                parent=self, branch_attrs=branch_attrs, root_name=root_name, deep=deep
            )
            for c in neighbors
        ]

        # identify node from which branch attributes should be transferred
        # 1. starting point (becomes root)
        if parent is None:
            other = None
        # 2. walk up (parent becomes child)
        elif parent.parent is self:
            other = parent
        # 3. walk down (retain the same order)
        else:
            other = self

        # create a new node and attach children to it
        result = self.__class__(children=children)

        # transfer attributes to the new node
        efc = self._exclude_from_copy
        for key in self.__dict__:
            if key not in efc:
                source = other if key in branch_attrs else self
                if source is not None and key in source.__dict__:
                    result.__dict__[key] = _copy(source.__dict__[key])

        # name the new root
        if root_name and parent is None and result.name is None:
            result.name = root_name

        return result

    def unrooted_deepcopy(self, parent=None):
        r"""Walk the tree unrooted-style and returns a new deepcopy.

        Parameters
        ----------
        parent : TreeNode or None
            Direction of walking (from parent to self). If specified, walking
            to the parent will be prohibited.

        Returns
        -------
        TreeNode
            A new copy of the tree rooted at the given node.

        Warnings
        --------
        ``unrooted_deepcopy`` is deprecated as of ``0.6.2``, as it generates a
        redundant copy of the tree. Use ``unrooted_copy`` instead.

        See Also
        --------
        copy
        unrooted_copy
        root_at

        Notes
        -----
        Perform a deepcopy of self and return a new copy of the tree as an
        unrooted copy. This is useful for defining a new root of the tree.

        This method calls ``unrooted_copy`` which is recursive.

        """
        msg = "Use unrooted_copy instead."
        _warn_deprecated(self.__class__.unrooted_deepcopy, "0.6.2", msg)

        root = self.root()
        root.assign_ids()

        new_tree = root.copy()
        new_tree.assign_ids()

        new_tree_self = new_tree.find_by_id(self.id)
        return new_tree_self.unrooted_copy(parent, deep=True)

    def unrooted_move(
        self,
        parent=None,
        branch_attrs={"length", "support"},
    ):
        r"""Walk the tree unrooted-style and rearrange it.

        .. versionadded:: 0.6.2

        Parameters
        ----------
        parent : TreeNode or None
            Direction of walking (from parent to self). If specified, walking
            to the parent will be prohibited.
        branch_attrs : set of str, optional
            Attributes of ``TreeNode`` objects that should be considered as
            branch attributes during the operation.

        See Also
        --------
        root_at
        unrooted_copy

        Notes
        -----
        This method recursively walks a tree from a given node in an unrooted
        style (i.e., directions of branches are not assumed). It rerranges the
        tree such that the given node becomes the root node and all other nodes
        are re-positioned accordingly, whereas the topology remains the same.

        This method manipulates the tree in place. There is no return value.
        The new tree should be referred to by the node where the operation
        started, as it has become the new root node.

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,(b,c)d)e,(f,g)h)i;"])
        >>> new_root = tree.find('d')
        >>> new_root.unrooted_move()
        >>> print(new_root)
        (b,c,(a,((f,g)h)i)e)d;
        <BLANKLINE>

        """
        # recursively add parent to children
        children = self.children
        if (old_parent := self.parent) is not None:
            children.append(old_parent)
            old_parent.unrooted_move(parent=self)

        # 1. starting point (becomes root)
        if parent is None:
            self.parent = None
            for attr in branch_attrs:
                setattr(self, attr, None)

        # 2. walk up (parent becomes child)
        else:
            for i, child in enumerate(children):
                if child is parent:
                    children.pop(i)
                    break
            self.parent = parent
            for attr in branch_attrs:
                setattr(self, attr, getattr(parent, attr, None))

    def count(self, tips=False):
        r"""Get the count of nodes in the tree.

        Parameters
        ----------
        tips : bool
            If ``True``, only return the count of tips.

        Returns
        -------
        int
            The number of nodes.

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,(b,c)d)e,(f,g)h)i;"])
        >>> print(tree.count())
        9
        >>> print(tree.count(tips=True))
        5

        """
        if tips:
            return len(list(self.tips()))
        else:
            return len(list(self.traverse(include_self=True)))

    def observed_node_counts(self, tip_counts):
        """Return counts of node observations from counts of tip observations.

        Parameters
        ----------
        tip_counts : dict of ints
            Counts of observations of tips. Keys correspond to tip names in
            ``self``, and counts are unsigned ints.

        Returns
        -------
        dict
            Counts of observations of nodes. Keys correspond to node names
            (internal nodes or tips), and counts are unsigned ints.

        Raises
        ------
        ValueError
            If a count less than one is observed.
        MissingNodeError
            If a count is provided for a tip not in the tree, or for an
            internal node.

        """
        result = defaultdict(int)
        for tip_name, count in tip_counts.items():
            if count < 1:
                raise ValueError("All tip counts must be greater than zero.")
            else:
                t = self.find(tip_name)
                if not t.is_tip():
                    raise MissingNodeError(
                        "Counts can only be for tips in the tree. %s is an "
                        "internal node." % t.name
                    )
                result[t] += count
                for internal_node in t.ancestors():
                    result[internal_node] += count
        return result

    def subtree(self, tip_list=None):
        r"""Make a copy of the subtree."""
        raise NotImplementedError()

    def subset(self):
        r"""Return set of names that descend from specified node.

        Get the set of `name` on tips that descend from this node.

        Returns
        -------
        frozenset
            The set of names at the tips of the clade that descends from self

        See Also
        --------
        subsets
        compare_subsets

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,(b,c)d)e,(f,g)h)i;"])
        >>> sorted(tree.subset())
        ['a', 'b', 'c', 'f', 'g']

        """
        return frozenset({i.name for i in self.tips()})

    def subsets(self):
        r"""Return all sets of names that come from self and its descendants.

        Compute all subsets of tip names over `self`, or, represent a tree as a
        set of nested sets.

        Returns
        -------
        frozenset
            A frozenset of frozensets of str

        See Also
        --------
        subset
        compare_subsets

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["(((a,b)c,(d,e)f)h)root;"])
        >>> subsets = tree.subsets()
        >>> len(subsets)
        3

        """
        sets = []
        for i in self.postorder(include_self=False):
            if not i.children:
                i.__leaf_set = frozenset([i.name])
            else:
                leaf_set = reduce(or_, [c.__leaf_set for c in i.children])
                if len(leaf_set) > 1:
                    sets.append(leaf_set)
                i.__leaf_set = leaf_set
        return frozenset(sets)

    def unroot(self, side=None):
        r"""Convert a rooted tree into unrooted.

        .. versionadded:: 0.6.2

        Parameters
        ----------
        side : int, optional
            Which basal node (i.e., children of root) will be elevated to root.
            Must be 0 or 1. If not provided, will elevate the first basal node
            that is not a tip.

        See Also
        --------
        root
        root_at

        Notes
        -----
        In scikit-bio, every tree has a root node. A tree is considered as
        "rooted" if its root node has exactly two children. In contrast, an
        "unrooted" tree may have three (the most common case), one, or more
        than three children attached to its root node. This method will not
        modify the tree if it is already unrooted.

        This method unroots a tree by trifucating its root. Specifically, it
        removes one of the two basal nodes of the tree (i.e., children of the
        root), transfers the name of the removed node to the root, and
        re-attaches the removed node's children to the root. Additionally, the
        removed node's branch length, if available, will be added to the other
        basal node's branch. The outcome appears as if the root is removed
        and the two basal nodes are directly connected.

        The choice of the basal node to be elevated affects the positioning of
        the resulting tree, but does not affect its topology from a
        phylogenetic perspective, as it is considered as unrooted.

        This method manipulates the tree in place. There is no return value.

        .. note:: In the case where the basal node has just one child, the
            resulting tree will still appear rooted as it has two basal nodes.
            To avoid this scenario, call ``prune`` to remove all one-child
            internal nodes.

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(['(((a,b)c,(d,e)f)g,(h,i)j)k;'])
        >>> print(tree.ascii_art())
                                      /-a
                            /c-------|
                           |          \-b
                  /g-------|
                 |         |          /-d
                 |          \f-------|
        -k-------|                    \-e
                 |
                 |          /-h
                  \j-------|
                            \-i

        >>> tree.unroot()
        >>> print(tree.ascii_art())
                            /-a
                  /c-------|
                 |          \-b
                 |
                 |          /-d
        -g-------|-f-------|
                 |          \-e
                 |
                 |          /-h
                  \j-------|
                            \-i

        """
        # return original tree if already unrooted
        root = self.root()
        if len(bases := root.children) != 2:
            return root

        # choose a basal node to elevate
        if side is None:
            side = 1 if (bases[0].is_tip() and not bases[1].is_tip()) else 0
        chosen, other = bases[side], bases[1 - side]

        # remove chosen node and re-attach its children to root
        root.invalidate_caches()
        chosen.parent = None
        for child in chosen.children:
            child.parent = root
        if side:
            root.children = [other] + chosen.children
        else:
            root.children = chosen.children + [other]

        # transfer basal node's name to root
        root.name = chosen.name
        # TODO: also transfer other custom node attributes

        # add branch length to the other basal node
        if (L := chosen.length) is not None:
            if other.length is not None:
                other.length += L
            else:
                other.length = L

    def _insert_above(self, above, branch_attrs=[]):
        """Insert a node into the branch connecting a node to its parent."""
        if above is False:
            return self
        node = self.__class__()
        if above is True:
            self.insert(node, None, branch_attrs)
        else:
            self.insert(node, above, branch_attrs)
        return node

    def root_at(
        self,
        node=None,
        above=False,
        reset=False,
        branch_attrs=["name"],
        root_name="root",
    ):
        r"""Reroot the tree at the provided node.

        This is useful for positioning a tree with an orientation that reflects
        knowledge of the true root location.

        Parameters
        ----------
        node : TreeNode or str, optional
            The node to root at. Can either be a node object or the name of the
            node. If not provided, will root at self. If a root node provided,
            will return the original tree.

            .. versionchanged:: 0.6.2

                Becomes optional.

        above : bool, float, or int, optional
            Whether and where to insert a new root node. If ``False``
            (default), the target node will serve as the root node. If
            ``True``, a new root node will be created and inserted at the
            midpoint of the branch connecting the target node and its parent.
            If a number, the new root will be inserted at this distance from
            the target node. The number ranges between 0 and branch length.

            .. versionadded:: 0.6.2

        reset : bool, optional
            Whether remove the original root of a rooted tree before performing
            the rerooting operation. Default is ``False``.

            .. versionadded:: 0.6.2

            .. note:: The default value will be set as ``True`` in 0.7.0.

        branch_attrs : iterable of str, optional
            Attributes of each node that should be considered as attributes of
            the branch connecting the node to its parent. This is important for
            the correct rerooting operation. "length" and "support" will be
            automatically included as they are always branch attributes.

            .. versionadded:: 0.6.2

            .. note:: ``name`` will be removed from the default in 0.7.0, as
                it is usually considered as an attribute of the node instead of
                the branch.

        root_name : str or None, optional
            Name for the root node, if it doesn't already have one.

            .. versionadded:: 0.6.2

            .. note:: The default value will be set as ``None`` in 0.7.0.

        Returns
        -------
        TreeNode
            A new copy of the tree rooted at the give node.

        Warnings
        --------
        The default behavior of ``root_at`` is subject to change in 0.7.0. The
        new default behavior can be achieved by specifying ``reset=True,
        branch_attrs=[], root_name=None``.

        See Also
        --------
        root_at_midpoint
        unrooted_copy
        unroot

        Notes
        -----
        The specified node will be come the root of the new tree.

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["(((a,b)c,(d,e)f)g,h)i;"])
        >>> print(tree.ascii_art())
                                      /-a
                            /c-------|
                           |          \-b
                  /g-------|
                 |         |          /-d
        -i-------|          \f-------|
                 |                    \-e
                 |
                  \-h

        Use the given node as the root node. This will typically create an
        unrooted tree (i.e., root node has three children).

        >>> t1 = tree.root_at("c", branch_attrs=[])
        >>> print(t1)
        (a,b,((d,e)f,(h)i)g)c;
        <BLANKLINE>
        >>> print(t1.ascii_art())
                  /-a
                 |
                 |--b
        -c-------|
                 |                    /-d
                 |          /f-------|
                  \g-------|          \-e
                           |
                            \i------- /-h

        Insert a new root node into the branch above the given node. This will
        create a rooted tree (i.e., root node has two children).

        >>> t2 = tree.root_at("c", above=True, branch_attrs=[])
        >>> print(t2)
        ((a,b)c,((d,e)f,(h)i)g)root;
        <BLANKLINE>
        >>> print(t2.ascii_art())
                            /-a
                  /c-------|
                 |          \-b
        -root----|
                 |                    /-d
                 |          /f-------|
                  \g-------|          \-e
                           |
                            \i------- /-h

        """
        # future warning
        if reset is False and branch_attrs == ["name"] and root_name == "root":
            func = self.__class__.root_at
            if not hasattr(func, "warned"):
                simplefilter("once", FutureWarning)
                warn(
                    "The default behavior of `root_at` is subject to change in 0.7.0. "
                    "The new default behavior can be achieved by specifying "
                    "`reset=True, branch_attrs=[], root_name=None`.",
                    FutureWarning,
                )
                func.warned = True

        tree = self.root()
        if node is None:
            node = self
        elif isinstance(node, str):
            node = tree.find(node)
        if node.is_root():
            return node.copy()

        if reset and len(tree.children) != 2:
            reset = False

        # copy the tree if it needs to be manipulated prior to walking
        if reset or above is not False:
            tree.assign_ids()
            new_tree = tree.copy()
            new_tree.assign_ids()
            node = new_tree.find_by_id(node.id)
            tree = new_tree

        # remove original root; we need to make sure the node itself is not the
        # basal node that gets removed
        if reset:
            side = None
            for i, base in enumerate(tree.children):
                if node is base:
                    side = 1 - i
                    break
            tree.unroot(side)

        # insert a new root node into the branch above
        node = node._insert_above(above, branch_attrs)

        branch_attrs = set(branch_attrs)
        branch_attrs.update(["length", "support"])
        return node.unrooted_copy(branch_attrs=branch_attrs, root_name=root_name)

    def root_at_midpoint(self, reset=False, branch_attrs=["name"], root_name="root"):
        r"""Reroot the tree at the midpoint of the two tips farthest apart.

        Parameters
        ----------
        reset : bool, optional
            Whether remove the original root of a rooted tree before performing
            the rerooting operation. Default is ``False``.

            .. versionadded:: 0.6.2

            .. note:: The default value will be set as ``True`` in 0.7.0.

        branch_attrs : iterable of str, optional
            Attributes of each node that should be considered as attributes of
            the branch connecting the node to its parent. This is important for
            the correct rerooting operation. "length" and "support" will be
            automatically included as they are always branch attributes.

            .. versionadded:: 0.6.2

            .. note:: ``name`` will be removed from the default in 0.7.0, as
                it is usually considered as an attribute of the node instead of
                the branch.

        root_name : str or None, optional
            Name for the new root node, if it doesn't have one.

            .. versionadded:: 0.6.2

            .. note:: The default value will be set as ``None`` in 0.7.0.

        Returns
        -------
        TreeNode
            A tree rooted at its midpoint.

        Raises
        ------
        TreeError
            If a tip ends up being the mid point.
        LengthError
            Midpoint rooting requires `length` and will raise (indirectly) if
            evaluated nodes don't have length.

        Warnings
        --------
        The default behavior of ``root_at_midpoint`` is subject to change in
        0.7.0. The new default behavior can be achieved by specifying
        ``reset=True, branch_attrs=[], root_name=None``.

        See Also
        --------
        root_at
        unrooted_copy

        Notes
        -----
        The midpoint rooting (MPR) method was originally described in [1]_.

        References
        ----------
        .. [1] Farris, J. S. (1972). Estimating phylogenetic trees from
           distance matrices. The American Naturalist, 106(951), 645-668.

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a:1,b:1)c:2,(d:3,e:4)f:5,g:1)h;"])
        >>> print(tree.ascii_art())
                            /-a
                  /c-------|
                 |          \-b
                 |
        -h-------|          /-d
                 |-f-------|
                 |          \-e
                 |
                  \-g

        >>> t = tree.root_at_midpoint(branch_attrs=[])
        >>> print(t)
        ((d:3.0,e:4.0)f:2.0,((a:1.0,b:1.0)c:2.0,g:1.0)h:3.0)root;
        <BLANKLINE>
        >>> print(t.ascii_art())
                            /-d
                  /f-------|
                 |          \-e
        -root----|
                 |                    /-a
                 |          /c-------|
                  \h-------|          \-b
                           |
                            \-g

        """
        # future warning
        if reset is False and branch_attrs == ["name"] and root_name == "root":
            func = self.__class__.root_at_midpoint
            if not hasattr(func, "warned"):
                simplefilter("once", FutureWarning)
                warn(
                    "The default behavior of `root_at_midpoint` is subject to change "
                    "in 0.7.0. The new default behavior can be achieved by specifying "
                    "`reset=True, branch_attrs=[], root_name=None`.",
                    FutureWarning,
                )
                func.warned = True

        tree = self.copy()
        if reset:
            tree.unroot()

        max_dist, tips = tree.get_max_distance()
        half_max_dist = max_dist / 2.0

        if max_dist == 0.0:
            return tree

        tip1 = tree.find(tips[0])
        tip2 = tree.find(tips[1])
        lca = tree.lowest_common_ancestor([tip1, tip2])

        if tip1.accumulate_to_ancestor(lca) > half_max_dist:
            climb_node = tip1
        else:
            climb_node = tip2

        dist_climbed = 0.0
        while dist_climbed + climb_node.length < half_max_dist:
            dist_climbed += climb_node.length
            climb_node = climb_node.parent

        # case 1: midpoint is at the climb node's parent
        # make the parent node as the new root
        if dist_climbed + climb_node.length == half_max_dist:
            new_root = climb_node.parent

        # case 2: midpoint is on the climb node's branch to its parent
        # insert a new root node into the branch
        else:
            new_root = tree.__class__()
            climb_node.insert(new_root, half_max_dist - dist_climbed)
            # TODO: Here, `branch_attrs` should be added to `insert`. However, this
            # will cause a backward-incompatible behavior. This change will be made
            # in version 0.7.0, along with the removal of `name` from the default of
            # `branch_attrs`.

        branch_attrs = set(branch_attrs)
        branch_attrs.update(["length", "support"])
        return new_root.unrooted_copy(branch_attrs=branch_attrs, root_name=root_name)

    def root_by_outgroup(
        self, outgroup, above=True, reset=True, branch_attrs=[], root_name=None
    ):
        r"""Reroot the tree with a given set of taxa as outgroup.

        .. versionadded:: 0.6.2

        Parameters
        ----------
        outgroup : iterable of str
            Taxon set to serve as outgroup. Must be a proper subset of taxa in
            the tree. The tree will be rooted at the lowest common ancestor
            (LCA) of the outgroup.
        above : bool, float, or int, optional
            Whether and where to insert a new root node. If ``False``, the
            LCA will serve as the root node. If ``True`` (default), a new root
            node will be created and inserted at the midpoint of the branch
            connecting the LCA and its parent (i.e., the midpoint between
            outgroup and ingroup). If a number between 0 and branch length, the
            new root will be inserted at this distance from the LCA.
        reset : bool, optional
            Whether remove the original root of a rooted tree before performing
            the rerooting operation. Default is ``True``.
        branch_attrs : iterable of str, optional
            Attributes of each node that should be considered as attributes of
            the branch connecting the node to its parent. This is important for
            the correct rerooting operation. "length" and "support" will be
            automatically included as they are always branch attributes.
        root_name : str or None, optional
            Name for the root node, if it doesn't already have one.

        Returns
        -------
        TreeNode
            A tree rooted by the outgroup.

        Raises
        ------
        TreeError
            Outgroup is not a proper subset of taxa in the tree.
        TreeError
            Outgroup is not monophyletic in the tree.

        Notes
        -----
        An outgroup is a subset of taxa that are usually distantly related from
        the remaining taxa (ingroup). The outgroup helps with locating the root
        of the ingroup, which are of interest in the study.

        This method reroots the tree at the lowest common ancestor (LCA) of the
        outgroup. By default, a new root will be placed at the midpoint between
        the LCA of outgroup and that of ingroup. But this behavior can be
        customized.

        This method requires the outgroup to be monophyletic, i.e., it forms a
        single clade in the tree. If the outgroup spans across the root of the
        tree, the method will reroot the tree within the ingroup such that the
        outgroup can form a clade in the rerooted tree, prior to rooting by
        outgroup.

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(['((((a,b),(c,d)),(e,f)),g);'])
        >>> print(tree.ascii_art())
                                                /-a
                                      /--------|
                                     |          \-b
                            /--------|
                           |         |          /-c
                           |          \--------|
                  /--------|                    \-d
                 |         |
                 |         |          /-e
        ---------|          \--------|
                 |                    \-f
                 |
                  \-g

        >>> rooted = tree.root_by_outgroup(['a', 'b'])
        >>> print(rooted.ascii_art())
                            /-a
                  /--------|
                 |          \-b
                 |
        ---------|                    /-c
                 |          /--------|
                 |         |          \-d
                  \--------|
                           |                    /-e
                           |          /--------|
                            \--------|          \-f
                                     |
                                      \-g

        >>> rooted = tree.root_by_outgroup(['e', 'f', 'g'])
        >>> print(rooted.ascii_art())
                                      /-e
                            /--------|
                  /--------|          \-f
                 |         |
                 |          \-g
        ---------|
                 |                    /-c
                 |          /--------|
                 |         |          \-d
                  \--------|
                           |          /-b
                            \--------|
                                      \-a

        """
        outgroup = set(outgroup)

        if not outgroup < self.subset():
            raise TreeError("Outgroup is not a proper subset of taxa in the tree.")

        # locate the lowest common ancestor (LCA) of outgroup in the tree
        lca = self.lca(outgroup)

        # if LCA is root (i.e., outgroup is split across basal clades), root
        # the tree at a tip within the ingroup and locate LCA again
        if lca is self:
            for tip in self.tips():
                if tip.name not in outgroup:
                    tree = self.root_at(tip, reset=reset, branch_attrs=branch_attrs)
                    break
            lca = tree.lca(outgroup)
        else:
            tree = self

        # test if outgroup is monophyletic
        if lca.count(tips=True) > len(outgroup):
            raise TreeError("Outgroup is not monophyletic in the tree.")

        # reroot the tree at LCA
        return tree.root_at(
            lca,
            above=above,
            reset=reset,
            branch_attrs=branch_attrs,
            root_name=root_name,
        )

    def is_tip(self):
        r"""Return `True` if the current node has no `children`.

        Returns
        -------
        bool
            `True` if the node is a tip

        See Also
        --------
        is_root
        has_children

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c);"])
        >>> print(tree.is_tip())
        False
        >>> print(tree.find('a').is_tip())
        True

        """
        return not self.children

    def is_root(self):
        r"""Return `True` if the current is a root, i.e. has no `parent`.

        Returns
        -------
        bool
            `True` if the node is the root

        See Also
        --------
        is_tip
        has_children

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c);"])
        >>> print(tree.is_root())
        True
        >>> print(tree.find('a').is_root())
        False

        """
        return self.parent is None

    def has_children(self):
        r"""Return `True` if the node has `children`.

        Returns
        -------
        bool
            `True` if the node has children.

        See Also
        --------
        is_tip
        is_root

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c);"])
        >>> print(tree.has_children())
        True
        >>> print(tree.find('a').has_children())
        False

        """
        return not self.is_tip()

    def traverse(self, self_before=True, self_after=False, include_self=True):
        r"""Return iterator over descendants.

        This is a depth-first traversal. Since the trees are not binary,
        preorder and postorder traversals are possible, but inorder traversals
        would depend on the data in the tree and are not handled here.

        Parameters
        ----------
        self_before : bool
            includes each node before its descendants if True
        self_after : bool
            includes each node after its descendants if True
        include_self : bool
            include the initial node if True


        `self_before` and `self_after` are independent. If neither is `True`,
        only terminal nodes will be returned.

        Note that if self is terminal, it will only be included once even if
        `self_before` and `self_after` are both `True`.

        Yields
        ------
        TreeNode
            Traversed node.

        See Also
        --------
        preorder
        postorder
        pre_and_postorder
        levelorder
        tips
        non_tips

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c);"])
        >>> for node in tree.traverse():
        ...     print(node.name)
        None
        c
        a
        b

        """
        if self_before:
            if self_after:
                return self.pre_and_postorder(include_self=include_self)
            else:
                return self.preorder(include_self=include_self)
        else:
            if self_after:
                return self.postorder(include_self=include_self)
            else:
                return self.tips(include_self=include_self)

    def preorder(self, include_self=True):
        r"""Perform preorder iteration over tree.

        Parameters
        ----------
        include_self : bool
            include the initial node if True

        Yields
        ------
        TreeNode
            Traversed node.

        See Also
        --------
        traverse
        postorder
        pre_and_postorder
        levelorder
        tips
        non_tips

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c);"])
        >>> for node in tree.preorder():
        ...     print(node.name)
        None
        c
        a
        b

        """
        stack = [self]
        while stack:
            curr = stack.pop()
            if include_self or (curr is not self):
                yield curr
            if curr.children:
                stack.extend(curr.children[::-1])

    def postorder(self, include_self=True):
        r"""Perform postorder iteration over tree.

        This is somewhat inelegant compared to saving the node and its index
        on the stack, but is 30% faster in the average case and 3x faster in
        the worst case (for a comb tree).

        Parameters
        ----------
        include_self : bool
            include the initial node if True

        Yields
        ------
        TreeNode
            Traversed node.

        See Also
        --------
        traverse
        preorder
        pre_and_postorder
        levelorder
        tips
        non_tips

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c);"])
        >>> for node in tree.postorder():
        ...     print(node.name)
        a
        b
        c
        None

        """
        child_index_stack = [0]
        curr = self
        curr_children = self.children
        curr_children_len = len(curr_children)
        while 1:
            curr_index = child_index_stack[-1]
            # if there are children left, process them
            if curr_index < curr_children_len:
                curr_child = curr_children[curr_index]
                # if the current child has children, go there
                if curr_child.children:
                    child_index_stack.append(0)
                    curr = curr_child
                    curr_children = curr.children
                    curr_children_len = len(curr_children)
                    curr_index = 0
                # otherwise, yield that child
                else:
                    yield curr_child
                    child_index_stack[-1] += 1
            # if there are no children left, return self, and move to
            # self's parent
            else:
                if include_self or (curr is not self):
                    yield curr
                if curr is self:
                    break
                curr = curr.parent
                curr_children = curr.children
                curr_children_len = len(curr_children)
                child_index_stack.pop()
                child_index_stack[-1] += 1

    def pre_and_postorder(self, include_self=True):
        r"""Perform iteration over tree, visiting node before and after.

        Parameters
        ----------
        include_self : bool
            include the initial node if True

        Yields
        ------
        TreeNode
            Traversed node.

        See Also
        --------
        traverse
        postorder
        preorder
        levelorder
        tips
        non_tips

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c);"])
        >>> for node in tree.pre_and_postorder():
        ...     print(node.name)
        None
        c
        a
        b
        c
        None

        """
        # handle simple case first
        if not self.children:
            if include_self:
                yield self
            return
        child_index_stack = [0]
        curr = self
        curr_children = self.children
        while 1:
            curr_index = child_index_stack[-1]
            if not curr_index:
                if include_self or (curr is not self):
                    yield curr
            # if there are children left, process them
            if curr_index < len(curr_children):
                curr_child = curr_children[curr_index]
                # if the current child has children, go there
                if curr_child.children:
                    child_index_stack.append(0)
                    curr = curr_child
                    curr_children = curr.children
                    curr_index = 0
                # otherwise, yield that child
                else:
                    yield curr_child
                    child_index_stack[-1] += 1
            # if there are no children left, return self, and move to
            # self's parent
            else:
                if include_self or (curr is not self):
                    yield curr
                if curr is self:
                    break
                curr = curr.parent
                curr_children = curr.children
                child_index_stack.pop()
                child_index_stack[-1] += 1

    def levelorder(self, include_self=True):
        r"""Perform levelorder iteration over tree.

        Parameters
        ----------
        include_self : bool
            include the initial node if True

        Yields
        ------
        TreeNode
            Traversed node.

        See Also
        --------
        traverse
        postorder
        preorder
        pre_and_postorder
        tips
        non_tips

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c,(d,e)f);"])
        >>> for node in tree.levelorder():
        ...     print(node.name)
        None
        c
        f
        a
        b
        d
        e

        """
        queue = [self]
        while queue:
            curr = queue.pop(0)
            if include_self or (curr is not self):
                yield curr
            if curr.children:
                queue.extend(curr.children)

    def tips(self, include_self=False):
        r"""Iterate over tips descended from `self`.

        Node order is consistent between calls and is ordered by a
        postorder traversal of the tree.

        Parameters
        ----------
        include_self : bool
            include the initial node if True

        Yields
        ------
        TreeNode
            Traversed node.

        See Also
        --------
        traverse
        postorder
        preorder
        pre_and_postorder
        levelorder
        non_tips

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c,(d,e)f);"])
        >>> for node in tree.tips():
        ...     print(node.name)
        a
        b
        d
        e

        """
        for n in self.postorder(include_self=include_self):
            if n.is_tip():
                yield n

    def non_tips(self, include_self=False):
        r"""Iterate over nontips descended from self.

        `include_self`, if `True` (default is False), will return the current
        node as part of non_tips if it is a non_tip. Node order is consistent
        between calls and is ordered by a postorder traversal of the tree.


        Parameters
        ----------
        include_self : bool
            include the initial node if True

        Yields
        ------
        TreeNode
            Traversed node.

        See Also
        --------
        traverse
        postorder
        preorder
        pre_and_postorder
        levelorder
        tips

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c,(d,e)f);"])
        >>> for node in tree.non_tips():
        ...     print(node.name)
        c
        f

        """
        for n in self.postorder(include_self):
            if not n.is_tip():
                yield n

    def invalidate_caches(self, attr=True):
        r"""Delete lookup and attribute caches.

        Parameters
        ----------
        attr : bool, optional
            If ``True``, invalidate attribute caches created by
            `TreeNode.cache_attr`.

        See Also
        --------
        create_caches
        cache_attr
        find

        """
        if not self.is_root():
            self.root().invalidate_caches()
        else:
            if hasattr(self, "_tip_cache"):
                delattr(self, "_tip_cache")
            if hasattr(self, "_non_tip_cache"):
                delattr(self, "_non_tip_cache")

            if hasattr(self, "_registered_caches") and attr:
                for node in self.traverse():
                    for cache in self._registered_caches:
                        if hasattr(node, cache):
                            delattr(node, cache)

    def create_caches(self):
        r"""Construct an internal lookup table to facilitate searching by name.

        Raises
        ------
        DuplicateNodeError
            The tip cache requires that names are unique (with the exception of
            names that are ``None``).

        See Also
        --------
        invalidate_caches
        cache_attr
        find

        Notes
        -----
        This method will not cache nodes whose name is ``None``. This method
        will raise ``DuplicateNodeError`` if a name conflict in the tips
        is discovered, but will not raise if on internal nodes. This is
        because, in practice, the tips of a tree are required to be unique
        while no such requirement holds for internal nodes.

        """
        if not self.is_root():
            self.root().create_caches()
        else:
            if hasattr(self, "_tip_cache") and hasattr(self, "_non_tip_cache"):
                return

            self.invalidate_caches(attr=False)

            tip_cache = {}
            non_tip_cache = defaultdict(list)

            for node in self.postorder():
                name = node.name

                if name is None:
                    continue

                if node.is_tip():
                    if name in tip_cache:
                        raise DuplicateNodeError(
                            f"Tip with name '{name}' already exists."
                        )

                    tip_cache[name] = node
                else:
                    non_tip_cache[name].append(node)

            self._tip_cache = tip_cache
            self._non_tip_cache = non_tip_cache

    def find_all(self, name):
        r"""Find all nodes that match `name`.

        The first call to `find_all` will cache all nodes in the tree on the
        assumption that additional calls to `find_all` will be made.

        Parameters
        ----------
        name : TreeNode or str
            The name or node to find. If `name` is `TreeNode` then all other
            nodes with the same name will be returned.

        Raises
        ------
        MissingNodeError
            Raises if the node to be searched for is not found

        Returns
        -------
        list of TreeNode
            The nodes found

        See Also
        --------
        find
        find_by_id
        find_by_func

        Examples
        --------
        >>> from skbio.tree import TreeNode
        >>> tree = TreeNode.read(["((a,b)c,(d,e)d,(f,g)c);"])
        >>> for node in tree.find_all('c'):
        ...     print(node.name, node.children[0].name, node.children[1].name)
        c a b
        c f g
        >>> for node in tree.find_all('d'):
        ...     print(node.name, str(node))
        d (d,e)d;
        <BLANKLINE>
        d d;
        <BLANKLINE>

        """
        root = self.root()

        # if what is being passed in looks like a node, just return it
        if isinstance(name, root.__class__):
            return [name]

        root.create_caches()

        tip = root._tip_cache.get(name, None)
        nodes = root._non_tip_cache.get(name, [])

        nodes.append(tip) if tip is not None else None

        if not nodes:
            raise MissingNodeError(f"Node '{name}' is not in self.")
        else:
            return nodes

    def find(self, name):
        r"""Find a node by name.

        Parameters
        ----------
        name : TreeNode or str
            The name of the node to find. If a ``TreeNode`` object is provided,
            then it is simply returned.

        Raises
        ------
        MissingNodeError
            Raises if the node to be searched for is not found.

        Returns
        -------
        TreeNode
            The found node.

        See Also
        --------
        find_all
        find_by_id
        find_by_func

        Notes
        -----
        The first call to ``find`` will cache all nodes in the tree on the
        assumption that additional calls to ``find`` will be made.

        ``find`` will first attempt to find the node in the tips. If it cannot
        find a corresponding tip, then it will search through the internal
        nodes of the tree. In practice, phylogenetic trees and other common
        trees in biology do not have unique internal node names. As a result,
        this find method will only return the first occurrence of an internal
        node encountered on a postorder traversal of the tree.

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c,(d,e)f);"])
        >>> print(tree.find('c').name)
        c

        """
        root = self.root()

        # if what is being passed in looks like a node, just return it
        if isinstance(name, root.__class__):
            return name

        root.create_caches()
        node = root._tip_cache.get(name, None)

        if node is None:
            node = root._non_tip_cache.get(name, [None])[0]

        if node is None:
            raise MissingNodeError("Node %s is not in self" % name)
        else:
            return node

    def find_by_id(self, node_id):
        r"""Find a node by `id`.

        This search method is based from the root.

        Parameters
        ----------
        node_id : int
            The `id` of a node in the tree

        Returns
        -------
        TreeNode
            The tree node with the matching id

        Notes
        -----
        This method does not cache id associations. A full traversal of the
        tree is performed to find a node by an id on every call.

        Raises
        ------
        MissingNodeError
            This method will raise if the `id` cannot be found

        See Also
        --------
        find
        find_all
        find_by_func

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c,(d,e)f);"])
        >>> print(tree.find_by_id(2).name)
        d

        """
        # if this method gets used frequently, then we should cache by ID
        # as well
        root = self.root()
        root.assign_ids()

        node = None
        for n in self.traverse(include_self=True):
            if n.id == node_id:
                node = n
                break

        if node is None:
            raise MissingNodeError("ID %d is not in self" % node_id)
        else:
            return node

    def find_by_func(self, func):
        r"""Find all nodes given a function.

        This search method is based on the current subtree, not the root.

        Parameters
        ----------
        func : a function
            A function that accepts a TreeNode and returns `True` or `False`,
            where `True` indicates the node is to be yielded

        Yields
        ------
        TreeNode
            Node found by `func`.

        See Also
        --------
        find
        find_all
        find_by_id

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c,(d,e)f);"])
        >>> func = lambda x: x.parent == tree.find('c')
        >>> [n.name for n in tree.find_by_func(func)]
        ['a', 'b']

        """
        for node in self.traverse(include_self=True):
            if func(node):
                yield node

    def ancestors(self):
        r"""Return all ancestors back to the root.

        This call will return all nodes in the path back to root, but does not
        include the node instance that the call was made from.

        Returns
        -------
        list of TreeNode
            The path, toward the root, from self

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c,(d,e)f)root;"])
        >>> [node.name for node in tree.find('a').ancestors()]
        ['c', 'root']

        """
        result = []
        curr = self
        while not curr.is_root():
            result.append(curr.parent)
            curr = curr.parent

        return result

    def root(self):
        r"""Return root of the tree which contains `self`.

        Returns
        -------
        TreeNode
            The root of the tree

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c,(d,e)f)root;"])
        >>> tip_a = tree.find('a')
        >>> root = tip_a.root()
        >>> root == tree
        True

        """
        curr = self
        while not curr.is_root():
            curr = curr.parent
        return curr

    def siblings(self):
        r"""Return all nodes that are `children` of `self` `parent`.

        This call excludes `self` from the list.

        Returns
        -------
        list of TreeNode
            The list of sibling nodes relative to self

        See Also
        --------
        neighbors

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c,(d,e,f)g)root;"])
        >>> tip_e = tree.find('e')
        >>> [n.name for n in tip_e.siblings()]
        ['d', 'f']

        """
        if self.is_root():
            return []

        result = self.parent.children[:]
        result.remove(self)

        return result

    def neighbors(self, ignore=None):
        r"""Return all nodes that are connected to self.

        This call does not include `self` in the result

        Parameters
        ----------
        ignore : TreeNode
            A node to ignore

        Returns
        -------
        list of TreeNode
            The list of all nodes that are connected to self

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c,(d,e)f)root;"])
        >>> node_c = tree.find('c')
        >>> [n.name for n in node_c.neighbors()]
        ['a', 'b', 'root']

        """
        nodes = [n for n in self.children + [self.parent] if n is not None]
        if ignore is None:
            return nodes
        else:
            return [n for n in nodes if n is not ignore]

    def lowest_common_ancestor(self, tipnames):
        r"""Find lowest common ancestor for a list of tips.

        Parameters
        ----------
        tipnames : iterable of TreeNode or str
            The nodes of interest

        Returns
        -------
        TreeNode
            The lowest common ancestor of the passed in nodes

        Raises
        ------
        ValueError
            If no tips could be found in the tree, or if not all tips were
            found.

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c,(d,e)f)root;"])
        >>> nodes = [tree.find('a'), tree.find('b')]
        >>> lca = tree.lowest_common_ancestor(nodes)
        >>> print(lca.name)
        c
        >>> nodes = [tree.find('a'), tree.find('e')]
        >>> lca = tree.lca(nodes)  # lca is an alias for convience
        >>> print(lca.name)
        root

        """
        if len(tipnames) == 1:
            return self.find(next(iter(tipnames)))

        tips = [self.find(name) for name in tipnames]

        if len(tips) == 0:
            raise ValueError("No tips found.")

        nodes_to_scrub = []

        for t in tips:
            if t.is_root():
                # has to be the LCA...
                return t

            prev = t
            curr = t.parent

            while curr and not hasattr(curr, "black"):
                setattr(curr, "black", [prev])
                nodes_to_scrub.append(curr)
                prev = curr
                curr = curr.parent

            # increase black count, multiple children lead to here
            if curr:
                curr.black.append(prev)

        curr = self
        while len(curr.black) == 1:
            curr = curr.black[0]

        # clean up tree
        for n in nodes_to_scrub:
            delattr(n, "black")

        return curr

    lca = lowest_common_ancestor  # for convenience

    @classonlymethod
    def from_taxonomy(cls, lineage_map):
        r"""Construct a tree from a taxonomy.

        Parameters
        ----------
        lineage_map : dict, iterable of tuples, or pd.DataFrame
            Mapping of taxon IDs to lineages (iterables of taxonomic units
            from high to low in ranking).

        Returns
        -------
        TreeNode
            The constructed taxonomy.

        See Also
        --------
        from_taxdump

        Examples
        --------
        >>> from skbio.tree import TreeNode
        >>> lineages = [
        ...     ('1', ['Bacteria', 'Firmicutes', 'Clostridia']),
        ...     ('2', ['Bacteria', 'Firmicutes', 'Bacilli']),
        ...     ('3', ['Bacteria', 'Bacteroidetes', 'Sphingobacteria']),
        ...     ('4', ['Archaea', 'Euryarchaeota', 'Thermoplasmata']),
        ...     ('5', ['Archaea', 'Euryarchaeota', 'Thermoplasmata']),
        ...     ('6', ['Archaea', 'Euryarchaeota', 'Halobacteria']),
        ...     ('7', ['Archaea', 'Euryarchaeota', 'Halobacteria']),
        ...     ('8', ['Bacteria', 'Bacteroidetes', 'Sphingobacteria']),
        ...     ('9', ['Bacteria', 'Bacteroidetes', 'Cytophagia'])]
        >>> tree = TreeNode.from_taxonomy(lineages)
        >>> print(tree.ascii_art())
                                      /Clostridia-1
                            /Firmicutes
                           |          \Bacilli- /-2
                  /Bacteria|
                 |         |                    /-3
                 |         |          /Sphingobacteria
                 |          \Bacteroidetes      \-8
                 |                   |
        ---------|                    \Cytophagia-9
                 |
                 |                              /-4
                 |                    /Thermoplasmata
                 |                   |          \-5
                  \Archaea- /Euryarchaeota
                                     |          /-6
                                      \Halobacteria
                                                \-7

        """
        root = cls(name=None)
        root._lookup = {}

        if isinstance(lineage_map, dict):
            lineage_map = lineage_map.items()
        elif isinstance(lineage_map, pd.DataFrame):
            lineage_map = ((idx, row.tolist()) for idx, row in lineage_map.iterrows())

        for id_, lineage in lineage_map:
            cur_node = root

            # for each name, see if we've seen it, if not, add that puppy on
            for name in lineage:
                if name in cur_node._lookup:
                    cur_node = cur_node._lookup[name]
                else:
                    new_node = cls(name=name)
                    new_node._lookup = {}
                    cur_node._lookup[name] = new_node
                    cur_node.append(new_node)
                    cur_node = new_node

            cur_node.append(cls(name=id_))

        # scrub the lookups
        for node in root.non_tips(include_self=True):
            del node._lookup

        return root

    def _balanced_distance_to_tip(self):
        """Return the distance to tip from this node.

        The distance to every tip from this node must be equal for this to
        return a correct result.

        Returns
        -------
        int
            The distance to tip of a length-balanced tree

        """
        node = self
        distance = 0
        while node.has_children():
            distance += node.children[0].length
            node = node.children[0]
        return distance

    @classonlymethod
    def from_linkage_matrix(cls, linkage_matrix, id_list):
        """Return tree from SciPy linkage matrix.

        Parameters
        ----------
        linkage_matrix : ndarray
            A SciPy linkage matrix as returned by
            `scipy.cluster.hierarchy.linkage`
        id_list : list
            The indices of the `id_list` will be used in the linkage_matrix

        Returns
        -------
        TreeNode
            An unrooted bifurcated tree

        See Also
        --------
        scipy.cluster.hierarchy.linkage

        """
        tip_width = len(id_list)
        cluster_count = len(linkage_matrix)
        lookup_len = cluster_count + tip_width
        node_lookup = np.empty(lookup_len, dtype=cls)

        for i, name in enumerate(id_list):
            node_lookup[i] = cls(name=name)

        for i in range(tip_width, lookup_len):
            node_lookup[i] = cls()

        newest_cluster_index = cluster_count + 1
        for link in linkage_matrix:
            child_a = node_lookup[int(link[0])]
            child_b = node_lookup[int(link[1])]

            path_length = link[2] / 2
            child_a.length = path_length - child_a._balanced_distance_to_tip()
            child_b.length = path_length - child_b._balanced_distance_to_tip()

            new_cluster = node_lookup[newest_cluster_index]
            new_cluster.append(child_a)
            new_cluster.append(child_b)

            newest_cluster_index += 1

        return node_lookup[-1]

    def to_taxonomy(self, allow_empty=False, filter_f=None):
        """Return a taxonomy representation of self.

        Parameters
        ----------
        allow_empty : bool, optional
            Allow gaps the taxonomy (e.g., internal nodes without names).
        filter_f : function, optional
            Specify a filtering function that returns True if the lineage is
            to be returned. This function must accept a ``TreeNode`` as its
            first parameter, and a ``list`` that represents the lineage as the
            second parameter.

        Yields
        ------
        tuple
            ``(tip, [lineage])`` where ``tip`` corresponds to a tip in the tree
            and ``[lineage]`` is the expanded names from root to tip. ``None``
            and empty strings are omitted from the lineage.

        Notes
        -----
        If ``allow_empty`` is ``True`` and the root node does not have a name,
        then that name will not be included. This is because it is common to
        have multiple domains represented in the taxonomy, which would result
        in a root node that does not have a name and does not make sense to
        represent in the output.

        Examples
        --------
        >>> from skbio.tree import TreeNode
        >>> lineages = {'1': ['Bacteria', 'Firmicutes', 'Clostridia'],
        ...             '2': ['Bacteria', 'Firmicutes', 'Bacilli'],
        ...             '3': ['Bacteria', 'Bacteroidetes', 'Sphingobacteria'],
        ...             '4': ['Archaea', 'Euryarchaeota', 'Thermoplasmata'],
        ...             '5': ['Archaea', 'Euryarchaeota', 'Thermoplasmata'],
        ...             '6': ['Archaea', 'Euryarchaeota', 'Halobacteria'],
        ...             '7': ['Archaea', 'Euryarchaeota', 'Halobacteria'],
        ...             '8': ['Bacteria', 'Bacteroidetes', 'Sphingobacteria'],
        ...             '9': ['Bacteria', 'Bacteroidetes', 'Cytophagia']}
        >>> tree = TreeNode.from_taxonomy(lineages.items())
        >>> lineages = sorted([(n.name, l) for n, l in tree.to_taxonomy()])
        >>> for name, lineage in lineages:
        ...     print(name, '; '.join(lineage))
        1 Bacteria; Firmicutes; Clostridia
        2 Bacteria; Firmicutes; Bacilli
        3 Bacteria; Bacteroidetes; Sphingobacteria
        4 Archaea; Euryarchaeota; Thermoplasmata
        5 Archaea; Euryarchaeota; Thermoplasmata
        6 Archaea; Euryarchaeota; Halobacteria
        7 Archaea; Euryarchaeota; Halobacteria
        8 Bacteria; Bacteroidetes; Sphingobacteria
        9 Bacteria; Bacteroidetes; Cytophagia

        """
        if filter_f is None:

            def filter_f(a, b):
                return True

        self.assign_ids()
        seen = set()
        lineage = []

        # visit internal nodes while traversing out to the tips, and on the
        # way back up
        for node in self.traverse(self_before=True, self_after=True):
            if node.is_tip():
                if filter_f(node, lineage):
                    yield (node, lineage[:])
            else:
                if allow_empty:
                    if node.is_root() and not node.name:
                        continue
                else:
                    if not node.name:
                        continue

                if node.id in seen:
                    lineage.pop(-1)
                else:
                    lineage.append(node.name)
                    seen.add(node.id)

    def to_array(self, attrs=None, nan_length_value=None):
        """Return an array representation of self.

        Parameters
        ----------
        attrs : list of tuple or None
            The attributes and types to return. The expected form is
            [(attribute_name, type)]. If `None`, then `name`, `length`, and
            `id` are returned.
        nan_length_value : float, optional
            If provided, replaces any `nan` in the branch length vector
            (i.e., ``result['length']``) with this value. `nan` branch lengths
            can arise from an edge not having a length (common for the root
            node parent edge), which can making summing problematic.

        Returns
        -------
        dict of array
            {id_index: {id: TreeNode},
             child_index: ((node_id, left_child_id, right_child_id)),
             attr_1: array(...),
             ...
             attr_N: array(...)}

        Notes
        -----
        Attribute arrays are in index order such that TreeNode.id can be used
        as a lookup into the array.

        Examples
        --------
        >>> from skbio import TreeNode
        >>> t = TreeNode.read(['(((a:1,b:2,c:3)x:4,(d:5)y:6)z:7);'])
        >>> res = t.to_array()
        >>> sorted(res.keys())
        ['child_index', 'id', 'id_index', 'length', 'name']
        >>> res['child_index'] # doctest: +ELLIPSIS
        array([[4, 0, 2],
               [5, 3, 3],
               [6, 4, 5],
               [7, 6, 6]]...
        >>> for k, v in res['id_index'].items():
        ...     print(k, v)
        ...
        0 a:1.0;
        <BLANKLINE>
        1 b:2.0;
        <BLANKLINE>
        2 c:3.0;
        <BLANKLINE>
        3 d:5.0;
        <BLANKLINE>
        4 (a:1.0,b:2.0,c:3.0)x:4.0;
        <BLANKLINE>
        5 (d:5.0)y:6.0;
        <BLANKLINE>
        6 ((a:1.0,b:2.0,c:3.0)x:4.0,(d:5.0)y:6.0)z:7.0;
        <BLANKLINE>
        7 (((a:1.0,b:2.0,c:3.0)x:4.0,(d:5.0)y:6.0)z:7.0);
        <BLANKLINE>
        >>> res['id']
        array([0, 1, 2, 3, 4, 5, 6, 7])
        >>> res['name']
        array(['a', 'b', 'c', 'd', 'x', 'y', 'z', None], dtype=object)

        """
        if attrs is None:
            attrs = [("name", object), ("length", float), ("id", int)]
        else:
            for attr, dtype in attrs:
                if not hasattr(self, attr):
                    raise AttributeError("Invalid attribute '%s'." % attr)

        id_index, child_index = self.index_tree()
        n = self.id + 1  # assign_ids starts at 0
        tmp = [np.zeros(n, dtype=dtype) for attr, dtype in attrs]

        for node in self.traverse(include_self=True):
            n_id = node.id
            for idx, (attr, dtype) in enumerate(attrs):
                tmp[idx][n_id] = getattr(node, attr)

        results = {"id_index": id_index, "child_index": child_index}
        results.update({attr: arr for (attr, dtype), arr in zip(attrs, tmp)})
        if nan_length_value is not None:
            length_v = results["length"]
            length_v[np.isnan(length_v)] = nan_length_value
        return results

    def _ascii_art(self, char1="-", show_internal=True, compact=False):
        LEN = 10
        PAD = " " * LEN
        PA = " " * (LEN - 1)
        namestr = self._node_label()
        if self.children:
            mids = []
            result = []
            for c in self.children:
                if c is self.children[0]:
                    char2 = "/"
                elif c is self.children[-1]:
                    char2 = "\\"
                else:
                    char2 = "-"
                (clines, mid) = c._ascii_art(char2, show_internal, compact)
                mids.append(mid + len(result))
                result.extend(clines)
                if not compact:
                    result.append("")
            if not compact:
                result.pop()
            (lo, hi, end) = (mids[0], mids[-1], len(result))
            prefixes = (
                [PAD] * (lo + 1) + [PA + "|"] * (hi - lo - 1) + [PAD] * (end - hi)
            )
            mid = int(np.trunc((lo + hi) / 2))
            prefixes[mid] = char1 + "-" * (LEN - 2) + prefixes[mid][-1]
            result = [p + L for (p, L) in zip(prefixes, result)]
            if show_internal:
                stem = result[mid]
                result[mid] = stem[0] + namestr + stem[len(namestr) + 1 :]
            return (result, mid)
        else:
            return ([char1 + "-" + namestr], 0)

    def ascii_art(self, show_internal=True, compact=False):
        r"""Return a string containing an ascii drawing of the tree.

        Note, this method calls a private recursive function and is not safe
        for large trees.

        Parameters
        ----------
        show_internal : bool
            includes internal edge names
        compact : bool
            use exactly one line per tip

        Returns
        -------
        str
            an ASCII formatted version of the tree

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b)c,(d,e)f)root;"])
        >>> print(tree.ascii_art())
                            /-a
                  /c-------|
                 |          \-b
        -root----|
                 |          /-d
                  \f-------|
                            \-e

        """
        (lines, mid) = self._ascii_art(show_internal=show_internal, compact=compact)
        return "\n".join(lines)

    def accumulate_to_ancestor(self, ancestor):
        r"""Return the sum of the distance between self and ancestor.

        Parameters
        ----------
        ancestor : TreeNode
            The node of the ancestor to accumulate distance too

        Returns
        -------
        float
            The sum of lengths between self and ancestor

        Raises
        ------
        NoParentError
            A NoParentError is raised if the ancestor is not an ancestor of
            self
        NoLengthError
            A NoLengthError is raised if one of the nodes between self and
            ancestor (including self) lacks a `length` attribute

        See Also
        --------
        distance

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a:1,b:2)c:3,(d:4,e:5)f:6)root;"])
        >>> root = tree
        >>> tree.find('a').accumulate_to_ancestor(root)
        4.0

        """
        accum = 0.0
        curr = self
        while curr is not ancestor:
            if curr.is_root():
                raise NoParentError("Provided ancestor is not in the path")

            if curr.length is None:
                raise NoLengthError(
                    "No length on node %s found." % curr.name or "unnamed"
                )

            accum += curr.length
            curr = curr.parent

        return accum

    def distance(self, other):
        """Return the distance between self and other.

        This method can be used to compute the distances between two tips,
        however, it is not optimized for computing pairwise tip distances.

        Parameters
        ----------
        other : TreeNode
            The node to compute a distance to

        Returns
        -------
        float
            The distance between two nodes

        Raises
        ------
        NoLengthError
            A NoLengthError will be raised if a node without `length` is
            encountered

        See Also
        --------
        tip_tip_distances
        accumulate_to_ancestor
        compare_tip_distances
        get_max_distance

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a:1,b:2)c:3,(d:4,e:5)f:6)root;"])
        >>> tip_a = tree.find('a')
        >>> tip_d = tree.find('d')
        >>> tip_a.distance(tip_d)
        14.0

        """
        if self is other:
            return 0.0

        self_ancestors = [self] + list(self.ancestors())
        other_ancestors = [other] + list(other.ancestors())

        if self in other_ancestors:
            return other.accumulate_to_ancestor(self)
        elif other in self_ancestors:
            return self.accumulate_to_ancestor(other)
        else:
            root = self.root()
            lca = root.lowest_common_ancestor([self, other])
            accum = self.accumulate_to_ancestor(lca)
            accum += other.accumulate_to_ancestor(lca)

            return accum

    def _set_max_distance(self):
        """Propagate tip distance information up the tree.

        This method was originally implemented by Julia Goodrich with the
        intent of being able to determine max tip to tip distances between
        nodes on large trees efficiently. The code has been modified to track
        the specific tips the distance is between
        """
        maxkey = itemgetter(0)

        for n in self.postorder():
            if n.is_tip():
                n.MaxDistTips = ((0.0, n), (0.0, n))
            else:
                if len(n.children) == 1:
                    raise TreeError("No support for single descedent nodes")
                else:
                    tip_info = [(max(c.MaxDistTips, key=maxkey), c) for c in n.children]

                    dists = [i[0][0] for i in tip_info]
                    best_idx = np.argsort(dists)[-2:]
                    (tip_a_d, tip_a), child_a = tip_info[best_idx[0]]
                    (tip_b_d, tip_b), child_b = tip_info[best_idx[1]]
                    tip_a_d += child_a.length or 0.0
                    tip_b_d += child_b.length or 0.0
                n.MaxDistTips = ((tip_a_d, tip_a), (tip_b_d, tip_b))

    def _get_max_distance_singledesc(self):
        """Return the max distance between any pair of tips.

        Also returns the tip names  that it is between as a tuple
        """
        distmtx = self.tip_tip_distances()
        idx_max = divmod(distmtx.data.argmax(), distmtx.shape[1])
        max_pair = (distmtx.ids[idx_max[0]], distmtx.ids[idx_max[1]])
        return distmtx[idx_max], max_pair

    def get_max_distance(self):
        """Return the max tip tip distance between any pair of tips.

        Returns
        -------
        float
            The distance between the two most distant tips in the tree
        tuple of TreeNode
            The two most distant tips in the tree

        Raises
        ------
        NoLengthError
            A NoLengthError will be thrown if a node without length is
            encountered

        See Also
        --------
        distance
        tip_tip_distances
        compare_tip_distances

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a:1,b:2)c:3,(d:4,e:5)f:6)root;"])
        >>> dist, tips = tree.get_max_distance()
        >>> dist
        16.0
        >>> [n.name for n in tips]
        ['b', 'e']

        """
        # _set_max_distance will throw a TreeError if a node with a single
        # child is encountered
        try:
            self._set_max_distance()
        except TreeError:  #
            return self._get_max_distance_singledesc()

        longest = 0.0
        tips = [None, None]
        for n in self.non_tips(include_self=True):
            tip_a, tip_b = n.MaxDistTips
            dist = tip_a[0] + tip_b[0]

            if dist > longest:
                longest = dist
                tips = [tip_a[1], tip_b[1]]

        # The MaxDistTips attribute causes problems during deep copy because it
        # contains references to other nodes. This patch removes the attribute.
        for n in self.traverse():
            del n.MaxDistTips

        return longest, tips

    def tip_tip_distances(self, endpoints=None):
        """Return distance matrix between pairs of tips, and a tip order.

        By default, all pairwise distances are calculated in the tree. If
        `endpoints` are specified, then only the distances between those tips
        are computed.

        Parameters
        ----------
        endpoints : list of TreeNode or str, or None
            A list of TreeNode objects or names of TreeNode objects

        Returns
        -------
        DistanceMatrix
            The distance matrix

        Raises
        ------
        ValueError
            If any of the specified `endpoints` are not tips

        See Also
        --------
        distance
        compare_tip_distances

        Notes
        -----
        If a node does not have an associated length, 0.0 will be used and a
        ``RepresentationWarning`` will be raised.

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a:1,b:2)c:3,(d:4,e:5)f:6)root;"])
        >>> mat = tree.tip_tip_distances()
        >>> print(mat)
        4x4 distance matrix
        IDs:
        'a', 'b', 'd', 'e'
        Data:
        [[  0.   3.  14.  15.]
         [  3.   0.  15.  16.]
         [ 14.  15.   0.   9.]
         [ 15.  16.   9.   0.]]

        """
        all_tips = list(self.tips())
        if endpoints is None:
            tip_order = all_tips
        else:
            tip_order = [self.find(n) for n in endpoints]
            for n in tip_order:
                if not n.is_tip():
                    raise ValueError("Node with name '%s' is not a tip." % n.name)

        # linearize all tips in postorder
        # .__start, .__stop compose the slice in tip_order.
        for i, node in enumerate(all_tips):
            node.__start, node.__stop = i, i + 1

        # the result map provides index in the result matrix
        result_map = {n.__start: i for i, n in enumerate(tip_order)}
        num_all_tips = len(all_tips)  # total number of tips
        num_tips = len(tip_order)  # total number of tips in result
        result = np.zeros((num_tips, num_tips), float)  # tip by tip matrix
        distances = np.zeros((num_all_tips), float)  # dist from tip to tip

        def update_result():
            # set tip_tip distance between tips of different child
            for child1, child2 in combinations(node.children, 2):
                for tip1 in range(child1.__start, child1.__stop):
                    if tip1 not in result_map:
                        continue
                    t1idx = result_map[tip1]
                    for tip2 in range(child2.__start, child2.__stop):
                        if tip2 not in result_map:
                            continue
                        t2idx = result_map[tip2]
                        result[t1idx, t2idx] = distances[tip1] + distances[tip2]

        for node in self.postorder():
            if not node.children:
                continue
            # subtree with solved child wedges
            # can possibly use np.zeros
            starts, stops = [], []  # to calc ._start and ._stop for curr node
            for child in node.children:
                length = child.length
                if length is None:
                    warn(
                        "`TreeNode.tip_tip_distances`: Node with name %r does "
                        "not have an associated length, so a length of 0.0 "
                        "will be used." % child.name,
                        RepresentationWarning,
                    )
                    length = 0.0
                distances[child.__start : child.__stop] += length

                starts.append(child.__start)
                stops.append(child.__stop)

            node.__start, node.__stop = min(starts), max(stops)

            if len(node.children) > 1:
                update_result()

        return DistanceMatrix(result + result.T, [n.name for n in tip_order])

    def compare_rfd(self, other, proportion=False):
        """Calculate the Robinson and Foulds symmetric difference.

        Parameters
        ----------
        other : TreeNode
            A tree to compare against
        proportion : bool
            Return a proportional difference

        Returns
        -------
        float
            The distance between the trees

        Notes
        -----
        Implementation based off of code by Julia Goodrich. The original
        description of the algorithm can be found in [1]_.

        Raises
        ------
        ValueError
            If the tip names between `self` and `other` are equal.

        See Also
        --------
        compare_subsets
        compare_tip_distances

        References
        ----------
        .. [1] Comparison of phylogenetic trees. Robinson and Foulds.
           Mathematical Biosciences. 1981. 53:131-141

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree1 = TreeNode.read(["((a,b),(c,d));"])
        >>> tree2 = TreeNode.read(["(((a,b),c),d);"])
        >>> tree1.compare_rfd(tree2)
        2.0

        """
        t1names = {n.name for n in self.tips()}
        t2names = {n.name for n in other.tips()}

        if t1names != t2names:
            if t1names < t2names:
                tree1 = self
                tree2 = other.shear(t1names)
            else:
                tree1 = self.shear(t2names)
                tree2 = other
        else:
            tree1 = self
            tree2 = other

        tree1_sets = tree1.subsets()
        tree2_sets = tree2.subsets()

        not_in_both = tree1_sets.symmetric_difference(tree2_sets)

        dist = float(len(not_in_both))

        if proportion:
            total_subsets = len(tree1_sets) + len(tree2_sets)
            dist /= total_subsets

        return dist

    def compare_subsets(self, other, exclude_absent_taxa=False):
        """Return fraction of overlapping subsets where self and other differ.

        Names present in only one of the two trees will count as mismatches,
        if you don't want this behavior, strip out the non-matching tips first.

        Parameters
        ----------
        other : TreeNode
            The tree to compare
        exclude_absent_taxa : bool
            Strip out names that don't occur in both trees

        Returns
        -------
        float
            The fraction of overlapping subsets that differ between the trees

        See Also
        --------
        compare_rfd
        compare_tip_distances
        subsets

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree1 = TreeNode.read(["((a,b),(c,d));"])
        >>> tree2 = TreeNode.read(["(((a,b),c),d);"])
        >>> tree1.compare_subsets(tree2)
        0.5

        """
        self_sets, other_sets = self.subsets(), other.subsets()

        if exclude_absent_taxa:
            in_both = self.subset() & other.subset()
            self_sets = (i & in_both for i in self_sets)
            self_sets = frozenset({i for i in self_sets if len(i) > 1})
            other_sets = (i & in_both for i in other_sets)
            other_sets = frozenset({i for i in other_sets if len(i) > 1})

        total_subsets = len(self_sets) + len(other_sets)
        intersection_length = len(self_sets & other_sets)

        if not total_subsets:  # no common subsets after filtering, so max dist
            return 1

        return 1 - (2 * intersection_length / float(total_subsets))

    def compare_tip_distances(
        self, other, sample=None, dist_f=distance_from_r, shuffle_f=np.random.shuffle
    ):
        """Compare self to other using tip-to-tip distance matrices.

        Value returned is `dist_f(m1, m2)` for the two matrices. Default is
        to use the Pearson correlation coefficient, with +1 giving a distance
        of 0 and -1 giving a distance of +1 (the maximum possible value).
        Depending on the application, you might instead want to use
        distance_from_r_squared, which counts correlations of both +1 and -1
        as identical (0 distance).

        Note: automatically strips out the names that don't match (this is
        necessary for this method because the distance between non-matching
        names and matching names is undefined in the tree where they don't
        match, and because we need to reorder the names in the two trees to
        match up the distance matrices).

        Parameters
        ----------
        other : TreeNode
            The tree to compare
        sample : int or None
            Randomly subsample the tips in common between the trees to
            compare. This is useful when comparing very large trees.
        dist_f : function
            The distance function used to compare two the tip-tip distance
            matrices
        shuffle_f : function
            The shuffling function used if `sample` is not None

        Returns
        -------
        float
            The distance between the trees

        Raises
        ------
        ValueError
            A ValueError is raised if there does not exist common tips
            between the trees

        See Also
        --------
        compare_subsets
        compare_rfd

        Examples
        --------
        >>> from skbio import TreeNode
        >>> # note, only three common taxa between the trees
        >>> tree1 = TreeNode.read(["((a:1,b:1):2,(c:0.5,X:0.7):3);"])
        >>> tree2 = TreeNode.read(["(((a:1,b:1,Y:1):2,c:3):1,Z:4);"])
        >>> dist = tree1.compare_tip_distances(tree2)
        >>> print("%.9f" % dist)
        0.000133446

        """
        self_names = {i.name: i for i in self.tips()}
        other_names = {i.name: i for i in other.tips()}
        common_names = frozenset(self_names) & frozenset(other_names)
        common_names = list(common_names)

        if not common_names:
            raise ValueError("No tip names in common between the two trees.")

        if len(common_names) <= 2:
            return 1  # the two trees must match by definition in this case

        if sample is not None:
            shuffle_f(common_names)
            common_names = common_names[:sample]

        self_nodes = [self_names[k] for k in common_names]
        other_nodes = [other_names[k] for k in common_names]

        self_matrix = self.tip_tip_distances(endpoints=self_nodes)
        other_matrix = other.tip_tip_distances(endpoints=other_nodes)

        return dist_f(self_matrix, other_matrix)

    def bifurcate(self, insert_length=None):
        r"""Reorder the tree into a bifurcating tree.

        All nodes that have more than two children will have additional
        intermediate nodes inserted to ensure that every node has only two
        children.

        Parameters
        ----------
        insert_length : int, optional
            The branch length assigned to all inserted nodes.

        See Also
        --------
        prune

        Notes
        -----
        Any nodes that have a single child can be collapsed using the
        prune method to create strictly bifurcating trees.

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b,g,h)c,(d,e)f)root;"])
        >>> print(tree.ascii_art())
                            /-a
                           |
                           |--b
                  /c-------|
                 |         |--g
                 |         |
        -root----|          \-h
                 |
                 |          /-d
                  \f-------|
                            \-e
        >>> tree.bifurcate()
        >>> print(tree.ascii_art())
                            /-h
                  /c-------|
                 |         |          /-g
                 |          \--------|
                 |                   |          /-a
        -root----|                    \--------|
                 |                              \-b
                 |
                 |          /-d
                  \f-------|
                            \-e

        """
        for n in self.traverse(include_self=True):
            if len(n.children) > 2:
                stack = n.children
                while len(stack) > 2:
                    ind = stack.pop()
                    intermediate = self.__class__()
                    intermediate.length = insert_length
                    intermediate.extend(stack)
                    n.append(intermediate)
                    for k in stack:
                        n.remove(k)
                    n.extend([ind, intermediate])

    def index_tree(self):
        """Index a tree for rapid lookups within a tree array.

        Indexes nodes in-place as `n._leaf_index`.

        Returns
        -------
        dict
            A mapping {node_id: TreeNode}
        np.array of ints
            This arrays describes the IDs of every internal node, and the ID
            range of the immediate descendents. The first column in the array
            corresponds to node_id. The second column is the left most
            descendent's ID. The third column is the right most descendent's
            ID.

        """
        self.assign_ids()

        id_index = {}
        child_index = []

        for n in self.postorder():
            for c in n.children:
                id_index[c.id] = c

                if c:
                    # c has children itself, so need to add to result
                    child_index.append((c.id, c.children[0].id, c.children[-1].id))

        # handle root, which should be t itself
        id_index[self.id] = self

        # only want to add to the child_index if self has children...
        if self.children:
            child_index.append((self.id, self.children[0].id, self.children[-1].id))
        child_index = np.asarray(child_index, dtype=np.int64)
        child_index = np.atleast_2d(child_index)

        return id_index, child_index

    def assign_ids(self):
        """Assign topologically stable unique ids to self.

        Following the call, all nodes in the tree will have their id
        attribute set.
        """
        curr_index = 0
        for n in self.postorder():
            for c in n.children:
                c.id = curr_index
                curr_index += 1

        self.id = curr_index

    def descending_branch_length(self, tip_subset=None):
        """Find total descending branch length from self or subset of self tips.

        Parameters
        ----------
        tip_subset : Iterable, or None
            If None, the total descending branch length for all tips in the
            tree will be returned. If a list of tips is provided then only the
            total descending branch length associated with those tips will be
            returned.

        Returns
        -------
        float
            The total descending branch length for the specified set of tips.

        Raises
        ------
        ValueError
            A ValueError is raised if the list of tips supplied to tip_subset
            contains internal nodes or non-tips.

        Notes
        -----
        This function replicates cogent's totalDescendingBranch Length method
        and extends that method to allow the calculation of total descending
        branch length of a subset of the tips if requested. The postorder
        guarantees that the function will always be able to add the descending
        branch length if the node is not a tip.

        Nodes with no length will have their length set to 0. The root length
        (if it exists) is ignored.

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tr = TreeNode.read(["(((A:.1,B:1.2)C:.6,(D:.9,E:.6)F:.9)G:2.4,"
        ...                     "(H:.4,I:.5)J:1.3)K;"])
        >>> tdbl = tr.descending_branch_length()
        >>> sdbl = tr.descending_branch_length(['A','E'])
        >>> print(round(tdbl, 1), round(sdbl, 1))
        8.9 2.2

        """
        self.assign_ids()
        if tip_subset is not None:
            all_tips = self.subset()
            if not set(tip_subset).issubset(all_tips):
                raise ValueError("tip_subset contains ids that aren't tip " "names.")

            lca = self.lowest_common_ancestor(tip_subset)
            ancestors = {}
            for tip in tip_subset:
                curr = self.find(tip)
                while curr is not lca:
                    ancestors[curr.id] = curr.length if curr.length is not None else 0.0
                    curr = curr.parent
            return sum(ancestors.values())

        else:
            return sum(
                n.length
                for n in self.postorder(include_self=False)
                if n.length is not None
            )

    def cache_attr(self, func, cache_attrname, cache_type=list):
        """Cache attributes on internal nodes of the tree.

        Parameters
        ----------
        func : function
            func will be provided the node currently being evaluated and must
            return a list of item (or items) to cache from that node or an
            empty list.
        cache_attrname : str
            Name of the attribute to decorate on containing the cached values
        cache_type : {set, frozenset, list}
            The type of the cache

        Notes
        -----
        This method is particularly useful if you need to frequently look up
        attributes that would normally require a traversal of the tree.

        WARNING: any cache created by this method will be invalidated if the
        topology of the tree changes (e.g., if `TreeNode.invalidate_caches` is
        called).

        Raises
        ------
        TypeError
            If an cache_type that is not a `set` or a `list` is specified.

        Examples
        --------
        Cache the tip names of the tree on its internal nodes

        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b,(c,d)e)f,(g,h)i)root;"])
        >>> f = lambda n: [n.name] if n.is_tip() else []
        >>> tree.cache_attr(f, 'tip_names')
        >>> for n in tree.traverse(include_self=True):
        ...     print("Node name: %s, cache: %r" % (n.name, n.tip_names))
        Node name: root, cache: ['a', 'b', 'c', 'd', 'g', 'h']
        Node name: f, cache: ['a', 'b', 'c', 'd']
        Node name: a, cache: ['a']
        Node name: b, cache: ['b']
        Node name: e, cache: ['c', 'd']
        Node name: c, cache: ['c']
        Node name: d, cache: ['d']
        Node name: i, cache: ['g', 'h']
        Node name: g, cache: ['g']
        Node name: h, cache: ['h']

        """
        if cache_type in (set, frozenset):

            def reduce_f(a, b):
                return a | b

        elif cache_type is list:

            def reduce_f(a, b):
                return a + b

        else:
            raise TypeError("Only list, set and frozenset are supported.")

        for node in self.postorder(include_self=True):
            if not hasattr(node, "_registered_caches"):
                node._registered_caches = set()
            node._registered_caches.add(cache_attrname)

            cached = [getattr(c, cache_attrname) for c in node.children]
            cached.append(cache_type(func(node)))
            setattr(node, cache_attrname, reduce(reduce_f, cached))

    def shuffle(self, k=None, names=None, shuffle_f=np.random.shuffle, n=1):
        """Yield trees with shuffled tip names.

        Parameters
        ----------
        k : int, optional
            The number of tips to shuffle. If k is not `None`, k tips are
            randomly selected, and only those names will be shuffled.
        names : list, optional
            The specific tip names to shuffle. k and names cannot be specified
            at the same time.
        shuffle_f : func
            Shuffle method, this function must accept a list and modify
            inplace.
        n : int, optional
            The number of iterations to perform. Value must be > 0 and `np.inf`
            can be specified for an infinite number of iterations.

        Notes
        -----
        Tip names are shuffled inplace. If neither `k` nor `names` are
        provided, all tips are shuffled.

        Yields
        ------
        TreeNode
            Tree with shuffled tip names.

        Raises
        ------
        ValueError
            If `k` is < 2
            If `n` is < 1
        ValueError
            If both `k` and `names` are specified
        MissingNodeError
            If `names` is specified but one of the names cannot be found

        Examples
        --------
        Alternate the names on two of the tips, 'a', and 'b', and do this 5
        times.

        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(["((a,b),(c,d));"])
        >>> rev = lambda items: items.reverse()
        >>> shuffler = tree.shuffle(names=['a', 'b'], shuffle_f=rev, n=5)
        >>> for shuffled_tree in shuffler:
        ...     print(shuffled_tree)
        ((b,a),(c,d));
        <BLANKLINE>
        ((a,b),(c,d));
        <BLANKLINE>
        ((b,a),(c,d));
        <BLANKLINE>
        ((a,b),(c,d));
        <BLANKLINE>
        ((b,a),(c,d));
        <BLANKLINE>

        """
        if k is not None and k < 2:
            raise ValueError("k must be None or >= 2")
        if k is not None and names is not None:
            raise ValueError("n and names cannot be specified at the sametime")
        if n < 1:
            raise ValueError("n must be > 0")

        self.assign_ids()

        if names is None:
            all_tips = list(self.tips())

            if n is None:
                n = len(all_tips)

            shuffle_f(all_tips)
            names = [tip.name for tip in all_tips[:k]]

        nodes = [self.find(name) for name in names]

        # Since the names are being shuffled, the association between ID and
        # name is no longer reliable
        self.invalidate_caches()

        counter = 0
        while counter < n:
            shuffle_f(names)
            for node, name in zip(nodes, names):
                node.name = name

            yield self
            counter += 1

    def _extract_support(self):
        """Extract the support value from a node label, if available.

        Returns
        -------
        tuple of
            int, float or None
                The support value extracted from the node label
            str or None
                The node label with the support value stripped

        """
        support, label = None, None
        if self.name:
            # separate support value from node name by the first colon
            left, _, right = self.name.partition(":")
            try:
                support = int(left)
            except ValueError:
                try:
                    support = float(left)
                except ValueError:
                    pass
            # strip support value from node name
            label = right or None if support is not None else self.name
        return support, label

    def _node_label(self):
        """Generate a node label.

        The label will be in the format of "support:name" if both exist,
        or "support" or "name" if either exists.

        Returns
        -------
        str
            Generated node label

        """
        lblst = []
        if self.support is not None:  # prevents support of NoneType
            lblst.append(str(self.support))
        if self.name:  # prevents name of NoneType
            lblst.append(self.name)
        return ":".join(lblst)

    def assign_supports(self):
        """Extract support values from internal node labels of a tree.

        Notes
        -----
        A "support value" measures the confidence or frequency of the incoming
        branch (the branch from parent to self) of an internal node in a tree.
        Roots and tips do not have support values. To extract a support value
        from a node label, this method reads from left and stops at the first
        ":" (if any), and attempts to convert it to a number.

        For examples: "(a,b)1.0", "(a,b)1.0:2.5", and "(a,b)'1.0:species_A'".
        In these cases the support values are all 1.0.

        For examples: "(a,b):1.0" and "(a,b)species_A". In these cases there
        are no support values.

        If a support value is successfully extracted, it will be stripped from
        the node label and assigned to the `support` property.

        IMPORTANT: mathematically, "support value" is a property of a branch,
        not a node, although they are usually attached to nodes in tree file
        formats [1]_.

        References
        ----------
        .. [1] Czech, Lucas, Jaime Huerta-Cepas, and Alexandros Stamatakis. "A
           Critical Review on the Use of Support Values in Tree Viewers and
           Bioinformatics Toolkits." Molecular biology and evolution 34.6
           (2017): 1535-1542.

        Examples
        --------
        >>> from skbio import TreeNode
        >>> newick = "((a,b)95,(c,d):1.1,(e,f)'80:speciesA':1.0);"
        >>> tree = TreeNode.read([newick])
        >>> tree.assign_supports()
        >>> tree.lca(['a', 'b']).support
        95
        >>> tree.lca(['c', 'd']).support is None
        True
        >>> tree.lca(['e', 'f']).support
        80
        >>> tree.lca(['e', 'f']).name
        'speciesA'

        """
        for node in self.traverse():
            if node.is_root() or node.is_tip():
                node.support = None
            else:
                node.support, node.name = node._extract_support()

    def unpack(self):
        """Unpack an internal node in place.

        Notes
        -----
        This method sequentially: 1) elongates child nodes by branch length
        of self (omit if there is no branch length), 2) removes self from
        parent node, and 3) grafts child nodes to parent node.

        Raises
        ------
        ValueError
            If input node is root or tip.

        See Also
        --------
        unpack_by_func
        prune

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(['((c:2.0,d:3.0)a:1.0,(e:2.0,f:1.0)b:2.0);'])
        >>> tree.find('b').unpack()
        >>> print(tree)
        ((c:2.0,d:3.0)a:1.0,e:4.0,f:3.0);
        <BLANKLINE>

        """
        if self.is_root():
            raise TreeError("Cannot unpack root.")
        if self.is_tip():
            raise TreeError("Cannot unpack tip.")
        parent = self.parent
        blen = self.length or 0.0
        for child in self.children:
            clen = child.length or 0.0
            child.length = clen + blen or None
        parent.remove(self)
        parent.extend(self.children)

    def unpack_by_func(self, func):
        """Unpack internal nodes of a tree that meet certain criteria.

        Parameters
        ----------
        func : function
            a function that accepts a TreeNode and returns `True` or `False`,
            where `True` indicates the node is to be unpacked

        See Also
        --------
        unpack
        prune

        Examples
        --------
        >>> from skbio import TreeNode
        >>> tree = TreeNode.read(['((c:2,d:3)a:1,(e:1,f:2)b:2);'])
        >>> tree.unpack_by_func(lambda x: x.length <= 1)
        >>> print(tree)
        ((e:1.0,f:2.0)b:2.0,c:3.0,d:4.0);
        <BLANKLINE>
        >>> tree = TreeNode.read(['(((a,b)85,(c,d)78)75,(e,(f,g)64)80);'])
        >>> tree.assign_supports()
        >>> tree.unpack_by_func(lambda x: x.support < 75)
        >>> print(tree)
        (((a,b)85,(c,d)78)75,(e,f,g)80);
        <BLANKLINE>

        """
        nodes_to_unpack = []
        for node in self.non_tips(include_self=False):
            if func(node):
                nodes_to_unpack.append(node)
        for node in nodes_to_unpack:
            node.unpack()

    @classonlymethod
    def from_taxdump(cls, nodes, names=None):
        r"""Construct a tree from the NCBI taxonomy database.

        Parameters
        ----------
        nodes : pd.DataFrame
            Taxon hierarchy
        names : pd.DataFrame or dict, optional
            Taxon names

        Returns
        -------
        TreeNode
            The constructed tree

        Notes
        -----
        ``nodes`` and ``names`` correspond to "nodes.dmp" and "names.dmp" of
        the NCBI taxonomy database. The should be read into data frames using
        ``skbio.io.read`` prior to this operation. Alternatively, ``names``
        may be provided as a dictionary. If ``names`` is omitted, taxonomy IDs
        be used as taxon names.

        Raises
        ------
        ValueError
            If there is no top-level node
        ValueError
            If there are more than one top-level node

        See Also
        --------
        from_taxonomy
        skbio.io.format.taxdump

        Examples
        --------
        >>> import pandas as pd
        >>> from skbio.tree import TreeNode
        >>> nodes = pd.DataFrame([
        ...             [1, 1, 'no rank'],
        ...             [2, 1, 'domain'],
        ...             [3, 1, 'domain'],
        ...             [4, 2, 'phylum'],
        ...             [5, 2, 'phylum']], columns=[
        ...     'tax_id', 'parent_tax_id', 'rank']).set_index('tax_id')
        >>> names = {1: 'root', 2: 'Bacteria', 3: 'Archaea',
        ...          4: 'Firmicutes', 5: 'Bacteroidetes'}
        >>> tree = TreeNode.from_taxdump(nodes, names)
        >>> print(tree.ascii_art())
                            /-Firmicutes
                  /Bacteria|
        -root----|          \-Bacteroidetes
                 |
                  \-Archaea

        """
        # identify top level of hierarchy
        tops = nodes[nodes["parent_tax_id"] == nodes.index]

        # validate root uniqueness
        n_top = tops.shape[0]
        if n_top == 0:
            raise ValueError("There is no top-level node.")
        elif n_top > 1:
            raise ValueError("There are more than one top-level node.")

        # get root taxid
        root_id = tops.index[0]

        # get parent-to-child(ren) map
        to_children = {
            p: g.index.tolist()
            for p, g in nodes[nodes.index != root_id].groupby("parent_tax_id")
        }

        # get rank map
        ranks = nodes["rank"].to_dict()

        # get taxon-to-name map
        # if not provided, use tax_id as name
        if names is None:
            names = {x: str(x) for x in nodes.index}

        # use "scientific name" as name
        elif isinstance(names, pd.DataFrame):
            names = names[names["name_class"] == "scientific name"][
                "name_txt"
            ].to_dict()

        # initiate tree
        tree = cls(names[root_id])
        tree.id = root_id
        tree.rank = ranks[root_id]

        # helper for extending tree
        def _extend_tree(node):
            self_id = node.id
            if self_id not in to_children:
                return
            children = []
            for id_ in to_children[self_id]:
                child = TreeNode(names[id_])
                child.id = id_
                child.rank = ranks[id_]
                _extend_tree(child)
                children.append(child)
            node.extend(children)

        # extend tree
        _extend_tree(tree)
        return tree
