TreeTime class documentation
TreeTime is the top-level wrapper class of the time tree inference package. In addition to inferring time trees, TreeTime can reroot your tree, resolve polytomies, mark tips that violate the molecular clock, or infer coalescent models. The core time tree inference is implemented in the class ClockTree.
TreeTime docstring and constructor
Main pipeline method
- TreeTime.run(root=None, infer_gtr=True, relaxed_clock=None, n_iqd=None, resolve_polytomies=True, max_iter=0, Tc=None, fixed_clock_rate=None, time_marginal='never', sequence_marginal=False, branch_length_mode='auto', vary_rate=False, use_covariation=False, tracelog_file=None, method_anc='probabilistic', assign_gamma=None, **kwargs)[source]
Run TreeTime reconstruction. Based on the input parameters, it divides the analysis into semi-independent jobs and conquers them one-by-one, gradually optimizing the tree given the temporal constarints and leaf node sequences.
- Parameters
root (str) –
Try to find better root position on a given tree. If string is passed, the root will be searched according to the specified method. If none, use tree as-is.
See
treetime.TreeTime.reroot()
for available rooting methods.infer_gtr (bool) – If True, infer GTR model
relaxed_clock (dict) – If not None, use autocorrelated molecular clock model. Specify the clock parameters as
{slack:<slack>, coupling:<coupling>}
dictionary.n_iqd (float) – If not None, filter tree nodes which do not obey the molecular clock for the particular tree. The nodes, which deviate more than
n_iqd
interquantile intervals from the molecular clock regression will be marked as ‘BAD’ and not used in the TreeTime analysisresolve_polytomies (bool) – If True, attempt to resolve multiple mergers
max_iter (int) – Maximum number of iterations to optimize the tree
Tc (float, str) –
If not None, use coalescent model to correct the branch lengths by introducing merger costs.
If Tc is float, it is interpreted as the coalescence time scale
If Tc is str, it should be one of (
opt
,const
,skyline
)fixed_clock_rate (float) – Fixed clock rate to be used. If None, infer clock rate from the molecular clock.
time_marginal (bool, str) – If False perform joint reconstruction of the divergence times, if True use marginal reconstruction of the divergence times, if ‘only_final’ (or ‘assign’) apply the marginal reconstruction only to the last optimization round, if “confidence-only” perform additional round using marginal reconstruction for calculation of confidence intervals but do not update times.
sequence_marginal (bool, optional) – use marginal reconstruction for ancestral sequences
branch_length_mode (str) –
Should be one of:
joint
,marginal
,input
.If ‘input’, rely on the branch lengths in the input tree and skip directly to the maximum-likelihood ancestral sequence reconstruction. Otherwise, perform preliminary sequence reconstruction using parsimony algorithm and do branch length optimization
vary_rate (bool or float, optional) – redo the time tree estimation for rates +/- one standard deviation. if a float is passed, it is interpreted as standard deviation, otherwise this standard deviation is estimated from the root-to-tip regression
use_covariation (bool, optional) – default False, if False, rate estimates will be performed using simple regression ignoring phylogenetic covaration between nodes. If vary_rate is True, use_covariation is true by default
method_anc (str, optional) – Which method should be used to reconstruct ancestral sequences. Supported values are “parsimony”, “fitch”, “probabilistic” and “ml”. Default is “probabilistic”
assign_gamma (callable, optional) – function to specify gamma (branch length scaling, local clock rate modifier) for each branch in tree, not compatible with a relaxed clock model
**kwargs – Keyword arguments needed by the downstream functions
- Returns
TreeTime error/succces code – return value depending on success or error
- Return type
str
Additional functionality
- TreeTime.resolve_polytomies(merge_compressed=False, resolution_threshold=0.05)[source]
Resolve the polytomies on the tree.
The function scans the tree, resolves polytomies if present, and re-optimizes the tree with new topology. Note that polytomies are only resolved if that would result in higher likelihood. Sometimes, stretching two or more branches that carry several mutations is less costly than an additional branch with zero mutations (long branches are not stiff, short branches are).
- Parameters
merge_compressed (bool) – If True, keep compressed branches as polytomies. If False, return a strictly binary tree.
- Returns
poly_found – The number of polytomies found
- Return type
int
- TreeTime.relaxed_clock(slack=None, coupling=None, **kwargs)[source]
Allow the mutation rate to vary on the tree (relaxed molecular clock). Changes of the mutation rates from one branch to another are penalized. In addition, deviation of the mutation rate from the mean rate is penalized.
- Parameters
slack (float) – Maximum change in substitution rate between parent and child nodes
coupling (float) – Maximum difference in substitution rates in sibling nodes
- TreeTime.clock_filter(reroot='least-squares', n_iqd=None, plot=False, fixed_clock_rate=None)[source]
Labels outlier branches that don’t seem to follow a molecular clock and excludes them from subsequent molecular clock estimation and the timetree propagation.
- Parameters
reroot (str) – Method to find the best root in the tree (see
treetime.TreeTime.reroot()
for options)n_iqd (float) –
Number of iqd intervals. The outlier nodes are those which do not fall into \(IQD\cdot n_iqd\) interval (\(IQD\) is the interval between 75th and 25th percentiles)
If None, the default (3) assumed
plot (bool) – If True, plot the results
- TreeTime.reroot(root='least-squares', force_positive=True, covariation=None, clock_rate=None)[source]
Find best root and re-root the tree to the new root
- Parameters
root (str) –
Which method should be used to find the best root. Available methods are:
best
, least-squares - minimize squared residual or likelihood of root-to-tip regressionmin_dev
- minimize variation of root-to-tip distanceoldest
- reroot on the oldest node<node_name>
- reroot to the node with name<node_name>
[<node_name1>, <node_name2>, ...]
- reroot to the MRCA of these nodes- force_positivebool
only consider positive rates when searching for the optimal root
- covariationbool
account for covariation in root-to-tip regression
- TreeTime.plot_root_to_tip(add_internal=False, label=True, ax=None)[source]
Plot root-to-tip regression
- Parameters
add_internal (bool) – If true, plot inte`rnal node positions
label (bool) – If true, label the plots
ax (matplotlib axes) – If not None, use the provided matplotlib axes to plot the results
- TreeTime.print_lh(joint=True)[source]
Print the total likelihood of the tree given the constrained leaves
- Parameters
joint (bool) – If true, print joint LH, else print marginal LH
- treetime.plot_vs_years(tt, step=None, ax=None, confidence=None, ticks=True, **kwargs)[source]
Converts branch length to years and plots the time tree on a time axis.
- Parameters
tt (TreeTime object) – A TreeTime instance after a time tree is inferred
step (int) – Width of shaded boxes indicating blocks of years. Will be inferred if not specified. To switch off drawing of boxes, set to 0
ax (matplotlib axes) – Axes to be used to plot, will create new axis if None
confidence (tuple, float) – Draw confidence intervals. This assumes that marginal time tree inference was run. Confidence intervals are either specified as an interval of the posterior distribution like (0.05, 0.95) or as the weight of the maximal posterior region , e.g. 0.9
**kwargs (dict) – Key word arguments that are passed down to Phylo.draw