dendropy.model.reconcile: Tree-Fitting (Gene/Species, Parasite/Host, etc.)

Classes and Methods for working with tree reconciliation, fitting, embedding, contained/containing etc.

class dendropy.model.reconcile.ContainingTree(containing_tree, contained_taxon_namespace, contained_to_containing_taxon_map, contained_trees=None, fit_containing_edge_lengths=True, collapse_empty_edges=True, ultrametricity_precision=False, ignore_root_deep_coalescences=True, **kwargs)

A “containing tree” is a (usually rooted) tree data structure within which other trees are “contained”. For example, species trees and their contained gene trees; host trees and their contained parasite trees; biogeographical “area” trees and their contained species or taxon trees.

__init__ converts self to ContainingTree class, embedding the trees given in the list, contained_trees.

Mandatory Arguments:

containing_tree
A Tree or Tree-like object that describes the topological constraints or conditions of the containing tree (e.g., species, host, or biogeographical area trees).
contained_taxon_namespace
A TaxonNamespace object that will be used to manage the taxa of the contained trees.
contained_to_containing_taxon_map
A TaxonNamespaceMapping object mapping Taxon objects in the contained TaxonNamespace to corresponding Taxon objects in the containing tree.

Optional Arguments:

contained_trees
An iterable container of Tree or Tree-like objects that will be contained into containing_tree; e.g. gene or parasite trees.
fit_containing_edge_lengths
If True [default], then the branch lengths of containing_tree will be adjusted to fit the contained tree as they are added. Otherwise, the containing tree edge lengths will not be changed.
collapse_empty_edges
If True [default], after edge lengths are adjusted, zero-length branches will be collapsed.
ultrametricity_precision
If False [default], then trees will not be checked for ultrametricity. Otherwise this is the threshold within which all node to tip distances for sister nodes must be equal.
ignore_root_deep_coalescences
If True [default], then deep coalescences in the root will not be counted.

Other Keyword Arguments: Will be passed to Tree().

build_edge_taxa_sets()

Rebuilds sets of containing and corresponding contained taxa at each edge.

clear()

Clears all contained trees and mapped edges.

clear_contained_edges()

Clears all contained mapped edges.

deep_coalescences()

Returns dictionary where the contained trees are keys, and the number of deep coalescences corresponding to the tree are values.

embed_contained_kingman(edge_pop_size_attr='pop_size', default_pop_size=1, label=None, rng=None, use_expected_tmrca=False)

Simulates, embeds, and returns a “censored” (Kingman) neutral coalescence tree conditional on self.

rng
Random number generator to use. If None, the default will be used.
edge_pop_size_attr
Name of attribute of self’s edges that specify the population size. If this attribute does not exist, then the population size is taken to be 1.

Note that all edge-associated taxon sets must be up-to-date (otherwise, build_edge_taxa_sets() should be called).

embed_tree(contained_tree)

Map edges of contained tree into containing tree (i.e., self).

fit_edge_lengths(contained_trees)

Recalculate node ages / edge lengths of containing tree to accomodate contained trees.

num_deep_coalescences()

Returns total number of deep coalescences of the contained trees.

rebuild(rebuild_taxa=True)

Recalculate edge taxa sets, node ages / edge lengths of containing tree, and embed edges of contained trees.

simulate_contained_kingman(edge_pop_size_attr='pop_size', default_pop_size=1, label=None, rng=None, use_expected_tmrca=False)

Simulates and returns a “censored” (Kingman) neutral coalescence tree conditional on self.

rng
Random number generator to use. If None, the default will be used.
edge_pop_size_attr
Name of attribute of self’s edges that specify the population size. If this attribute does not exist, then the population size is taken to be 1.

Each coalescent tree terminal node will have a container_tree_node attribute added that references the node on the container tree from which it was sampled.

Note that all edge-associated taxon sets must be up-to-date (otherwise, build_edge_taxa_sets() should be called), and that the tree is not added to the set of contained trees. For the latter, call embed_contained_kingman.

write_as_mesquite(out, **kwargs)

For debugging purposes, write out a Mesquite-format file.

dendropy.model.reconcile.monophyletic_partition_discordance(tree, taxon_namespace_partition)

Returns the number of deep coalescences on tree tree that would result if the taxa in tax_sets formed K mutually-exclusive monophyletic groups, where K = len(tax_sets) taxon_namespace_partition == TaxonNamespacePartition

dendropy.model.reconcile.reconciliation_discordance(gene_tree, species_tree)

Given two trees (with splits encoded), this returns the number of gene duplications implied by the gene tree reconciled on the species tree, based on the algorithm described here:

Goodman, M. J. Czelnusiniak, G. W. Moore, A. E. Romero-Herrera, and G. Matsuda. 1979. Fitting the gene lineage into its species lineage, a parsimony strategy illustrated by cladograms constructed from globin sequences. Syst. Zool. 19: 99-113.

Maddison, W. P. 1997. Gene trees in species trees. Syst. Biol. 46: 523-536.

This function requires that the gene tree and species tree have the same leaf set. Note that for correct results,

  1. trees must be rooted (i.e., is_rooted = True)
  2. split masks must have been added as rooted (i.e., when encode_splits was called, is_rooted must have been set to True)