ivy.chars package

Submodules

ivy.chars.bayesian_models module

ivy.chars.bayesian_models.create_mk_model(tree, chars, Qtype, pi)[source]

Create model objects to be passed to pymc.MCMC

Creates Qparams and likelihood function

ivy.chars.bayesian_models.create_multi_mk_model(tree, chars, Qtype, pi, nregime=2)[source]

Create an mk model with multiple regimes to be sampled from with MCMC.

Regime number is fixed and the location of the regime shift is allowed to change

ivy.chars.bayesian_models.fit_mk_bayes(tree, chars, Qtype, pi, *kwargs)[source]

Fit an mk model to a given tree and list of characters. Return posterior distributions of Q parameters and MAP estimate of Q matrix

Parameters:
  • tree (Node) – Root node of a tree. All branch lengths must be greater than 0 (except root)
  • chars (list) – List of character states corresponding to leaf nodes in preoder sequence. Character states must be in the form of 0,1,2,...
  • pi (str) – Either “Equal”, “Equilibrium”, or “Fitzjohn”. How to weight values at root node. Defaults to “Equal” Method “Fitzjohn” is not thouroughly tested, use with caution
  • Qtype –

    Either a string specifying how to esimate values for Q or a numpy array of a pre-specified Q matrix.

    Valid strings for Q:

    “Equal”: All rates equal “Sym”: Forward and reverse rates equal “ARD”: All rates different
Keyword Arguments:
 
  • iters (float) – Number of iterations in MCMC. Defaults to 2000
  • burn (float) – Burnin to discard. Defaults to 200
  • thin (float) – Thinning parameter. Defaults to 1
Returns:

The pymc MCMC object and the pymc MAP object
Return type:

tuple
ivy.chars.bayesian_models.get_indices(list_of_lists)[source]

Get indices given list of nodes in regimes

ivy.chars.bayesian_models.hrm_bayesian(tree, chars, Qtype, nregime, pi=u'Fitzjohn', constraint=u'Rate')[source]

Create a hidden rates model for pymc to be sampled from.

Parameters:
  • tree (Node) – Root node of a tree. All branch lengths must be greater than 0 (except root)
  • chars (dict) –

    Dict mapping character states to tip labels. Character states should be coded 0,1,2...

    Can also be a list with tip states in preorder sequence

  • pi (str) – Either “Equal”, “Equilibrium”, or “Fitzjohn”. How to weight values at root node. Defaults to “Equal” Method “Fitzjohn” is not thouroughly tested, use with caution
  • Qtype –

    Either a string specifying how to esimate values for Q or a numpy array of a pre-specified Q matrix.

    “Simple”: Symmetric rates within observed states and between
    rates.
    “STD”: State Transitions Different. Transitions between states
    within the same rate class are asymetrical
  • nregime (int) – Number of hidden states. nstates = 0 is equivalent to a vanilla Mk model
  • constraint (str) –

    Contraints to apply to the parameters of the Q matrix. Can be one of the following:

    “Rate”: The fastest rate in the fastest regime must be faster than
    the fastest rate in the slowest regime
    “Symmetry”: For two-regime models only. The two regimes must
    have different symmetry (a>b in regime 1, b>a in regime 2)

    “None”: No contraints
ivy.chars.bayesian_models.mk_results(mcmc_obj)[source]

Create summary graphs and statistics for an mk model

Parameters:mcmc_obj – A pymc MCMC object that has been sampled
Returns:Trace plots, histograms, and summary statistics
Return type:dict
ivy.chars.bayesian_models.nshifts(node, inds, n=0)[source]

The number of regime shifts, given indices of regimes

ivy.chars.bayesian_models.valid_indices(nobschar, nregime, x, y)[source]

Return only the valid indices for the given subarray

ivy.chars.catpars module

ivy.chars.catpars.ancstates(tree, chardata, stepmatrix)[source]

Return parsimony ancestral states

Parameters:
  • tree (Node) – Root of tree
  • chardata (Dict) – Dict of tip labels mapping to character states
  • stepmatrix (np.array) – Step matrix. Create default step matrix with default_costmatrix().
Returns:

Internal nodes mapping to reconstructed states.
Return type:

dict
ivy.chars.catpars.binary_deltran(tree, chardata, stepmatrix)[source]
ivy.chars.catpars.default_costmatrix(numstates, dtype=<type 'int'>)[source]

a square array with zeroes along the diagonal, ones elsewhere

ivy.chars.catpars.downpass(node, states, stepmatrix, chardata, node2dpv=None)[source]
ivy.chars.catpars.minstates(v)[source]

return the indices of v that equal the minimum

ivy.chars.catpars.uppass(node, states, stepmatrix, node2dpv, node2upm={}, node2ancstates=None)[source]

ivy.chars.discrete module

ivy.chars.discrete.NoTO(tree, chars)[source]

Number of Tips Per Origin

See: Bromham et al. 2016

ivy.chars.discrete.is_monotypic(node, chardict)[source]
ivy.chars.discrete.monotypic_clade_size(tree, chars)[source]

Count diversity contained within subclades having the same character state.

ivy.chars.discrete.nodeLikelihood(node)[source]

Take node “node” and calculate its likelihood given its children’s likelihoods, branch lengths, and p-matrix. :param node: A node to calculate the likelihood for :type node: Node

Returns:The likelihood of the node given the data
Return type:float
ivy.chars.discrete.sse_get_lambda(params, i, j, k, nstate)[source]
ivy.chars.discrete.sse_get_mu(params, i, nstate)[source]
ivy.chars.discrete.sse_get_qij(params, i, j, nstate)[source]
ivy.chars.discrete.tip_age_rank_sum(tree, chars)[source]

Calculate tip age rank sums of two traits and return test statistic and p-value

See: Bromham et al. 2016

ivy.chars.evolve module

Functions for evolving traits and trees.

ivy.chars.evolve.brownian(root, sigma=1.0, init=0.0, values={})[source]

Recursively evolve a trait by Brownian motion up from the node root.

  • sigma: standard deviation of the normal random variate after one unit of branch length
  • init: initial value

Returns: values - a dictionary mapping nodes to evolved values

ivy.chars.evolve.test_brownian()[source]

Evolve a trait up an example tree of primates:.

((((Homo:0.21,Pongo:0.21)N1:0.28,Macaca:0.49)N2:0.13, Ateles:0.62)N3:0.38,Galago:1.00)root;

Returns: (root, data) - the root node and evolved data.

ivy.chars.hrm module

ivy.chars.hrm.AIC(l, k)[source]

Akaike information criterion

ivy.chars.hrm.any_equal(ar)[source]

Given a 3-D matrix, test if any of the sub-matrices are equal. For use in evaluating symmetric constraint in hrm likelihood function

ivy.chars.hrm.cluster_models(tree, chars, Q, nregime, pi=u'Equal', findmin=True)[source]

Given an MLE Q, return candidate models with more parsimonious parameter set by merging similar parameters

ivy.chars.hrm.create_hrm_ar(tree, chars, nregime, findmin=True)[source]

Create arrays to be used in hrm likelihood function

ivy.chars.hrm.create_likelihood_function_hrm_mk_MLE(tree, chars, nregime, Qtype, pi=u'Equal', findmin=True, constraints=u'Rate', orderedRegimes=True)[source]

Create a function that takes values for Q and returns likelihood.

Specify the Q to be ER, Sym, or ARD

Returned function to be passed into scipy.optimize

Parameters:
  • tree (Node) – Root node of a tree. All branch lengths must be greater than 0 (except root)
  • chars (dict) –

    Dict mapping character states to tip labels. Character states should be coded 0,1,2...

    Can also be a list with tip states in preorder sequence

  • Qtype (str) – ARD only
  • pi (str) – Either “Equal”, “Equilibrium”, or “Fitzjohn”. How to weight values at root node.
  • min (bool) – Whether the function is to be minimized (False means it will be maximized)
  • orderedRegimes (bool) – Whether or not regimes are ordered such that a transition from 1 -> 3 is not possible (a branch in state 1 must pass through 2 to get to 3)

Notes

Constraints:
The maximum rate within each rate category must be ordered ( fastest rate in slower regime must be slower than fastest rate in fastest regime)
Returns:
Function accepting a list of parameters and returning
log-likelihood. To be optmimized with scipy.optimize.minimize
Return type:function
ivy.chars.hrm.extract_Qparams(Q, nregime)[source]
ivy.chars.hrm.extract_regimes(Q, nobschar, nregime)[source]

Given a HRM Q matrix,extract the within-regime transitions

ivy.chars.hrm.fill_Q_layout(regimetype, Qparams)[source]

Fill a single regime matrix based on the regime type and the given parameters.

ivy.chars.hrm.fill_Q_matrix(nobschar, nregime, wrparams, brparams, Qtype=u'ARD', out=None, orderedRegimes=True)[source]

Fill a Q matrix with nchar*nregime rows and cols with values from Qparams

Parameters:
  • nchar (int) – number of observed characters
  • nregime (int) – number of hidden rates per character
  • wrparams (list) – List of unique Q values for within-regime transitions, in order as they appear in columnwise iteration
  • brparams (list) – list of unique Q values for between-regime transition, in order as they appear in columnwise iteration
  • Qtype (str) – Either “ARD” or “Simple”. What type of Q matrix to fill.
  • out (np.array) – Optional numpy array to put output into
  • orderedRegimes (bool) – Whether or not regimes are ordered such that a transition from 1 -> 3 is not possible (a branch in state 1 must pass through 2 to get to 3)
Returns:

Q-matrix with values filled in. Check to make sure values

have been filled in properly
Return type:

array
ivy.chars.hrm.fill_distinct_regime_Q(regimetypes, Qparams, nregime, nobschar, Q=None, Q_layout=None, br_variable=False)[source]
ivy.chars.hrm.fill_model_Q(mod, Qparams, Q)[source]

Fill Q matrix given model and parameters.

Parameters:
  • mod (list) –

    List of tuples. Code for the mode used. The first nregime tuples correspond to the within-regime transition rates for each regime. The last tuple corresponds to the between-regime transition rates. Ex. [(1, 2), (3, 4), (5, 6, 7, 8)] corresponds to a matrix of:

    [-,1,5,-] [2,-,-,6] [7,-,-,3] [-,8,4,-]
  • Qparams (list) – List of floats corresponding to the values for slow, medium, fast, etc. rates. The first is always 0
  • Q (np.array) – Pre-allocated Q matrix
ivy.chars.hrm.fit_hrm(tree, chars, nregime, pi=u'Equal', constraints=u'Rate', Qtype=u'ARD', orderedRegimes=True, startingvals=None)[source]

Fit a hidden-rates mk model to a given tree and list of characters, and number of regumes. Return fitted ARD Q matrix and calculated likelihood.

Parameters:
  • tree (Node) – Root node of a tree. All branch lengths must be greater than 0 (except root)
  • chars (dict) –

    Dict mapping character states to tip labels. Character states should be coded 0,1,2...

    Can also be a list with tip states in preorder sequence

  • nregime (int) – Number of hidden rates per character
  • pi (str) – Either “Equal”, “Equilibrium”, or “Fitzjohn”. How to weight values at root node. Defaults to “Equal” Method “Fitzjohn” is not thouroughly tested, use with caution
  • orderedRegimes (bool) – Whether or not to constrain regime transitions to only occur between adjacent regimes
  • startingvals (list) – Values to initialize optimization at. Optional
Returns:

Tuple of fitted Q matrix (a np array) and log-likelihood value
Return type:

tuple
ivy.chars.hrm.fit_hrm_distinct_regimes(tree, chars, nregime, nparams, pi=u'Equal', br_variable=False, parallel=False, ncores=2, out_file=None)[source]

Fit hrm with distinct regime types given number of regimes and number of parameters

BINARY CHARACTERS ONLY

Parameters:
  • tree (Node) – Root node of tree
  • chars (dict) –

    Dict mapping character states to tip labels. Character states should be coded 0,1,2...

    Can also be a list with tip states in preorder sequence

  • nregime (int) – Number of regimes to test in model
  • nparams (int) – Number of unique parameters available for regimes
  • pi (str) – Root prior
  • br_variable (bool) – Whether or not between-regime rates are allowed to vary
  • parallel (bool) – Whether to run in parallel
  • ncores (int) – If parallel is True, number of cores to run on
  • out_file (str) – Optional output file name. If given, results will written into the current directory with this string as the first part of the filename
ivy.chars.hrm.fit_hrm_mkARD(tree, chars, nregime, pi=u'Equal', constraints=u'Rate', orderedRegimes=True, startingvals=None)[source]

Fit a hidden-rates mk model to a given tree and list of characters, and number of regumes. Return fitted ARD Q matrix and calculated likelihood.

Parameters:
  • tree (Node) – Root node of a tree. All branch lengths must be
  • chars (dict) –

    Dict mapping character states to tip labels. Character states should be coded 0,1,2...

    Can also be a list with tip states in preorder sequence

  • nregime (int) – Number of hidden rates per character
  • pi (str) – Either “Equal”, “Equilibrium”, or “Fitzjohn”. How to weight values at root node. Defaults to “Equal” Method “Fitzjohn” is not thouroughly tested, use with caution
  • orderedRegimes (bool) – Whether or not to constrain regime transitions to only occur between adjacent regimes
Returns:

Tuple of fitted Q matrix (a np array) and log-likelihood value
Return type:

tuple
ivy.chars.hrm.fit_hrm_mkSimple(tree, chars, nregime, pi=u'Equal', orderedRegimes=True, startingvals=None)[source]

Fit a hidden-rates mk model to a given tree and list of characters, and number of regumes. Return fitted ARD Q matrix and calculated likelihood.

Parameters:
  • tree (Node) – Root node of a tree. All branch lengths must be greater than 0 (except root)
  • chars (dict) –

    Dict mapping character states to tip labels. Character states should be coded 0,1,2...

    Can also be a list with tip states in preorder sequence

  • nregime (int) – Number of hidden rates per character
  • pi (str) – Either “Equal”, “Equilibrium”, or “Fitzjohn”. How to weight values at root node. Defaults to “Equal” Method “Fitzjohn” is not thouroughly tested, use with caution
Returns:

Tuple of fitted Q matrix (a np array) and log-likelihood value
Return type:

tuple
ivy.chars.hrm.fit_hrm_model(tree, chars, nregime, mod, pi=u'Equal', findmin=True, initialvals=None)[source]

Fit parameters for a single model specified by the user

ivy.chars.hrm.fit_hrm_model_likelihood(tree, chars, nregime, mod, pi=u'Equal', findmin=True)[source]

Likelihood function for fitting user-specified model

ivy.chars.hrm.fit_hrm_qidx(tree, chars, nregime, qidx, pi=u'Equal', orderedRegimes=True, startingvals=None)[source]

Fit a hidden-rates mk model to a given tree and list of characters, and number of regumes. Return fitted ARD Q matrix and calculated likelihood.

Parameters:
  • tree (Node) – Root node of a tree. All branch lengths must be
  • chars (dict) –

    Dict mapping character states to tip labels. Character states should be coded 0,1,2...

    Can also be a list with tip states in preorder sequence

  • nregime (int) – Number of hidden rates per character
  • pi (str) – Either “Equal”, “Equilibrium”, or “Fitzjohn”. How to weight values at root node. Defaults to “Equal” Method “Fitzjohn” is not thouroughly tested, use with caution
  • orderedRegimes (bool) – Whether or not to constrain regime transitions to only occur between adjacent regimes
Returns:

Dict of fitted Q matrix (a np array) and log-likelihood value
Return type:

dict
ivy.chars.hrm.fit_hrm_randstartingpoints(tree, chars, nregime, nstart=10, pi=u'Equal', constraints=u'Rate', Qtype=u'ARD', orderedRegimes=True)[source]

Fit a hidden-rates mk model to a given tree and list of characters, and number of regumes. Return fitted ARD Q matrix and calculated likelihood. Run nstart times, each with different randomized starting values, and return the highest likelihood

Parameters:
  • tree (Node) – Root node of a tree. All branch lengths must be greater than 0 (except root)
  • chars (dict) –

    Dict mapping character states to tip labels. Character states should be coded 0,1,2...

    Can also be a list with tip states in preorder sequence

  • nregime (int) – Number of hidden rates per character
  • nstart (int) – Number of times to run likelihood fit
  • pi (str) – Either “Equal”, “Equilibrium”, or “Fitzjohn”. How to weight values at root node. Defaults to “Equal” Method “Fitzjohn” is not thouroughly tested, use with caution
  • orderedRegimes (bool) – Whether or not to constrain regime transitions to only occur between adjacent regimes
Returns:

Tuple of fitted Q matrix (a np array) and log-likelihood value
Return type:

tuple
ivy.chars.hrm.format_mod(mod, nregime, nobschar)[source]

Convert flattened model to formatted model for use in fill_model_Q

ivy.chars.hrm.hrm_disctinct_regimes_likelihoodfunc(tree, chars, regimetypes, pi=u'Equal', findmin=True, br_variable=False, ar=None)[source]
ivy.chars.hrm.hrm_mk(tree, chars, Q, nregime, pi=u'Equal', returnPi=False, ar=None)[source]

Return log-likelihood of hidden-rates-markov mk model as described in Beaulieu et al. 2013

Parameters:
  • tree (Node) – Root node of a tree. All branch lengths must be greater than 0 (except root)
  • chars (dict) –

    Dict mapping character states to tip labels. Character states should be coded 0,1,2...

    Can also be a list with tip states in preorder sequence

  • Q (np.array) – Instantaneous rate matrix
  • p (np.array) – Optional pre-allocated p matrix
  • pi (str or np.array) –

    Option to weight the root node by given values. Either a string containing the method or an array of weights. Weights should be given in order. Accepted methods of weighting root:

    Equal: flat prior Equilibrium: Prior equal to stationary distribution
    of Q matrix
    Fitzjohn: Root states weighted by how well they
    explain the data at the tips.
  • returnPi (bool) – Whether or not to return the final values of root node weighting
  • ar (dict) – Dict of pre-allocated arrays to improve speed by avoiding creating and destroying new arrays
Returns:

Log-likelihood of model
Return type:

float

Examples

from ivy.examples import primates,primate_data Q = np.array([[-0.11,0.1 ,0.01 ,0 ],

[0.1 ,-0.11,0.0 ,0.01 ], [0.01 ,0.0 ,-0.51,0.5 ], [0.0 ,0.01 ,0.5 ,-0.51]])

hrm_mk(primates,primate_data,Q,nregime=2)

ivy.chars.hrm.make_identity_regime(regimePair)[source]

Given one regime pair, reduce it to its identity (eg (2,2),(3,3) becomes (1,1),(2,2)

ivy.chars.hrm.make_regime_type_combos(nregime, nparams)[source]

Create regime combinations for a binary character

ivy.chars.hrm.optimize_regime(comb, tree=None, chars=None, nregime=None, nparams=None, pi=None, br_variable=None, out_file=None, ar=None)[source]

Top-level function for optimizing hrm distinct-regime model in parallel.

ivy.chars.hrm.pairwise_merge(tree, chars, Q, nregime, pi=u'Equal')[source]

Given an MLE Q, merge similar parameters pairwise to find more parsimonious models

Parameters:
  • tree (Node) – Root node of a tree. All branch lengths must be greater than 0 (except root)
  • chars (dict) –

    Dict mapping character states to tip labels. Character states should be coded 0,1,2...

    Can also be a list with tip states in preorder sequence

  • Q (np.array) – Instantaneous rate matrix
  • nregime (int) – Number of regimes in the model
  • pi (str or np.array) –

    Option to weight the root node by given values. Either a string containing the method or an array of weights. Weights should be given in order. Accepted methods of weighting root:

    Equal: flat prior Equilibrium: Prior equal to stationary distribution
    of Q matrix
    Fitzjohn: Root states weighted by how well they
    explain the data at the tips.
ivy.chars.hrm.remove_redundant_regimes(regimePairs)[source]

Remove redundant regime pairs (eg (1,1),(2,2) and (2,2),(3,3) are redundant and only one should be kept)

ivy.chars.hrm_bayesian module

class ivy.chars.hrm_bayesian.QmatMetropolis(stochastic, nparam, nchar, nregime)[source]

Bases: pymc.StepMethods.Metropolis

Custom step algorithm for selecting a new model

propose()[source]
reject()[source]
class ivy.chars.hrm_bayesian.QmatMetropolis_mk(stochastic, all_mods)[source]

Bases: pymc.StepMethods.Metropolis

Custom step algorithm for selecting a new model

propose()[source]
reject()[source]
ivy.chars.hrm_bayesian.ShiftedGamma(shape, shift=1, name=u'ShiftedGamma')[source]

Gamma distribution that has been shifted by some constant.

ivy.chars.hrm_bayesian.fill_model_Q(mod, Qparams, Q)[source]

Create Q matrix given model and parameters.

Parameters:
  • mod (list) –

    List of tuples. Code for the mode used. The first nregime tuples correspond to the within-regime transition rates for each regime. The last tuple corresponds to the between-regime transition rates. Ex. [(1, 2), (3, 4), (5, 6, 7, 8)] corresponds to a matrix of:

    [-,1,5,-] [2,-,-,6] [7,-,-,3] [-,8,4,-]
  • Qparams (list) – List of floats corresponding to the values for slow, medium, fast, etc. rates. The first is always 0
  • Q (np.array) – Pre-allocated Q matrix
ivy.chars.hrm_bayesian.fill_model_mk(mod, Qparams, Q, mask)[source]

Fill cells of Q with Qparams based on mod

ivy.chars.hrm_bayesian.hrm_allmodels_bayes(tree, chars, nregime, nparam, modseed, pi=u'Equal', dbname=None)[source]

Use an MCMC chain to fit a hrm model with a limited number of parameters.

The chain will step through different models, placing a prior on simpler models

Parameters:
  • tree (Node) – Root node of a tree. All branch lengths must be greater than 0 (except root)
  • chars (dict) –

    Dict mapping character states to tip labels. Character states should be coded 0,1,2...

    Can also be a list with tip states in preorder sequence

  • nregime (int) – Number of regimes
  • nparam (int) – Number of unique parameters to allow in a model
  • modseed (tuple) – Starting model for the MCMC chain. A tuple of ints. Must be a valid model or adjacent to at least one valid model.
  • pi (str) – Either “Equal”, “Equilibrium”, or “Fitzjohn”. How to weight values at root node. Defaults to “Equal” Method “Fitzjohn” is not thouroughly tested, use with caution
ivy.chars.hrm_bayesian.is_valid_model(mod, nparam, nchar, nregime, mod_order)[source]

Check if a given model is valid

ivy.chars.hrm_bayesian.make_qmat_stoch(nobschar, nregime, nparam, mod_order_list, modseed, name=u'qmat_stoch')[source]

Make a stochastic to use for a model-changing step

ivy.chars.hrm_bayesian.make_qmat_stoch_mk(all_mods, name=u'qmat_stoch')[source]

Stochastic for use in model-changing step

ivy.chars.hrm_bayesian.mk_allmodels_bayes(tree, chars, nparam, pi=u'Equal', dbname=None)[source]

Fit an mk model with nparam parameters distributed about the Q matrix.

ivy.chars.hrm_bayesian.new_hrm_model(mod, nparam, nchar, nregime, mod_order)[source]

Return new, valid model by changing one parameter of current model

ivy.chars.hrm_bayesian.number_connected(mod, nchar, nregime, n_wr, n_br)[source]

Test how many regimes are connected in some way

ivy.chars.hrm_bayesian.unique_models(nchar, nregime, nparam)[source]

Create a list of all possible models for a Q matrix with nchar distinct character states and nregime distinct regimes, filling all cells with nparam + 0.

Parameters:
  • nchar (int) – Number of observed characters
  • nregime (int) – Number of regimes in model
  • nparam (int) – Number of unique parameters in model (not including 0)
Returns:

List of all unique models. Each item contains a tuple of the

within-regime parameters and a tuple of the between-regime parameters.
Return type:

list

ivy.chars.mk module

ivy.chars.mk.create_likelihood_function_mk(tree, chars, Qtype, pi=u'Equal', findmin=True)[source]

Create a function that takes values for Q and returns likelihood.

Specify the Q to be ER, Sym, or ARD

Returned function to be passed into scipy.optimize

Parameters:
  • tree (Node) – Root node of a tree. All branch lengths must be greater than 0 (except root)
  • chars (dict) –

    Dict mapping character states to tip labels. Character states should be coded 0,1,2...

    Can also be a list with tip states in preorder sequence

  • Qtype (str) – What type of Q matrix to use. Either ER (equal rates), Sym (symmetric rates), or ARD (All rates different).
  • pi (str) – Either “Equal”, “Equilibrium”, or “Fitzjohn”. How to weight values at root node.
  • min (bool) – Whether the function is to be minimized (False means it will be maximized)
Returns:

Function accepting a list of parameters and returning

log-likelihood. To be optmimized with NLOPT
Return type:

function
ivy.chars.mk.create_mk_ar(tree, chars, findmin=True)[source]

Create preallocated arrays. For use in mk function

Nodelist = edgelist of nodes in postorder sequence

ivy.chars.mk.fitMkARD(tree, chars, pi=u'Equal')[source]

Estimate parameters of an all-rates-different Q matrix Return log-likelihood of mk equation using fitted Q

Multi-parameter model: all rates different

Parameters:
  • tree (Node) – Root node of a tree. All branch lengths must be greater than 0 (except root)
  • chars (dict) –

    Dict mapping character states to tip labels. Character states should be coded 0,1,2...

    Can also be a list with tip states in preorder sequence

  • pi (str) – Either “Equal” or “Fitzjohn”. How to weight values at root node. Defaults to “Equal” Method “Fitzjohn” is currently untested
Returns:

Fitted parameter, log-likelihood, and dictionary of weightings

at the root.
Return type:

tuple
ivy.chars.mk.fitMkER(tree, chars, pi=u'Equal')[source]

Estimate parameter of an equal-rate Q matrix Return log-likelihood of mk equation using fitted Q

One-parameter model: alpha = beta

Parameters:
  • tree (Node) – Root node of a tree. All branch lengths must be greater than 0 (except root)
  • chars (dict) –

    Dict mapping character states to tip labels. Character states should be coded 0,1,2...

    Can also be a list with tip states in preorder sequence

  • pi (str) – Either “Equal” or “Fitzjohn”. How to weight values at root node. Defaults to “Equal”
Returns:

Fitted parameter, log-likelihood, and dictionary of weightings

at the root.
Return type:

tuple
ivy.chars.mk.fitMkSym(tree, chars, pi=u'Equal')[source]

Estimate parameter of a symmetrical-rate Q matrix Return log-likelihood of mk equation using fitted Q

Multi-parameter model: forward = reverse

Parameters:
  • tree (Node) – Root node of a tree. All branch lengths must be greater than 0 (except root)
  • chars (dict) –

    Dict mapping character states to tip labels. Character states should be coded 0,1,2...

    Can also be a list with tip states in preorder sequence

  • pi (str) – Either “Equal” or “Fitzjohn”. How to weight values at root node. Defaults to “Equal” Method “Fitzjohn” is currently untested
Returns:

Fitted parameter, log-likelihood, and dictionary of weightings

at the root.
Return type:

tuple
ivy.chars.mk.fit_Mk(tree, chars, Q=u'Equal', pi=u'Equal')[source]

Fit an mk model to a given tree and list of characters. Return fitted Q matrix and calculated likelihood.

Parameters:
  • tree (Node) – Root node of a tree. All branch lengths must be greater than 0 (except root)
  • chars (dict) –

    Dict mapping character states to tip labels. Character states should be coded 0,1,2...

    Can also be a list with tip states in preorder sequence

  • pi – Either “Equal”, “Equilibrium”, or “Fitzjohn”. How to weight values at root node. Defaults to “Equal” Method “Fitzjohn” is not thouroughly tested, use with caution
Returns:

Log-likelihood, fitted Q matrix, root prior, root likelihood
Return type:

dict

Examples

from ivy.examples import primates, primate_data primate_eq = fit_Mk(primates,primate_data,Q=”Equal”) primate_ARD = fit_Mk(primates, primate_data,Q=”ARD”)

ivy.chars.mk.mk(tree, chars, Q, p=None, pi=u'Equal', returnPi=False, ar=None)[source]

Fit mk model and return likelihood for the root node of a tree given a list of characters and a Q matrix

Convert tree and character data into a form that can be input into mk, which fits an mk model.

Note: internal calculations are log-transformed to avoid numeric underflow

Parameters:
  • tree (Node) – Root node of a tree. All branch lengths must be greater than 0 (except root)
  • chars (dict) –

    Dict mapping character states to tip labels. Character states should be coded 0,1,2...

    Can also be a list with tip states in preorder sequence

  • Q (np.array) – Instantaneous rate matrix
  • p (np.array) – Optional pre-allocated p matrix
  • pi (str or np.array) –

    Option to weight the root node by given values. Either a string containing the method or an array of weights. Weights should be given in order.

    Accepted methods of weighting root:

    Equal: flat prior Equilibrium: Prior equal to stationary distribution

    of Q matrix
    Fitzjohn: Root states weighted by how well they
    explain the data at the tips.
  • returnPi (bool) – Whether or not to return the final values of root node weighting
  • ar (dict) – Dict of pre-allocated arrays to improve speed by avoiding creating and destroying new arrays
Returns:

log-likelihood of model
Return type:

float

Examples

from ivy.examples import primates, primate_data Q = np.array([[-0.1,0.1],[0.1,-0.1]]) mk(primates,primate_data,Q)

ivy.chars.mk.qsd(Q)[source]

Calculate stationary distribution of Q, assuming each state has the same diversification rate.

Parameters:Q (np.array) – Instantaneous rate matrix
Returns:Stationary distribution of pi
Return type:(np.array)

Eqn from Maddison et al 2007

Referenced from Phytools (Revell 2013)

ivy.chars.mk_mr module

class ivy.chars.mk_mr.ModelOrderMetropolis(stochastic)[source]

Bases: pymc.StepMethods.Metropolis

Custom step method for model order

competence()[source]
propose()[source]
reject()[source]
class ivy.chars.mk_mr.SwitchpointMetropolis(stochastic, tree, seg_map, stepsize=0.05, seglen=0.02)[source]

Bases: pymc.StepMethods.Metropolis

Custom step algorithm for selecting a new switchpoint

competence()[source]
global_step()[source]
local_step()[source]
propose()[source]
reject()[source]
ivy.chars.mk_mr.create_likelihood_function_multimk(tree, chars, Qtype, nregime, pi=u'Equal', findmin=True)[source]
ivy.chars.mk_mr.create_likelihood_function_multimk_mods(tree, chars, mods, pi=u'Equal', findmin=True, orderedparams=True)[source]

Create a likelihood function for testing the parameter values of user- specified models

ivy.chars.mk_mr.create_mkmr_ar(tree, chars, nregime, findmin=True)[source]

Create preallocated arrays. For use in mk function

Nodelist = edgelist of nodes in postorder sequence

ivy.chars.mk_mr.create_mkmr_mb_ar(tree, chars, nregime, findmin=True)[source]

Preallocated arrays for midbranch mk_mr

Create preallocated arrays. For use in mk_mr midbranch function

Nodelist = edgelist of nodes in postorder sequence

ivy.chars.mk_mr.fill_model_mr_Q(Qparams, mods, Q)[source]

Fill a Q-matrix with Qparams corresponding to the model.

Individual Q matrices are filled rowwise

Function alters Q in place

ivy.chars.mk_mr.inds_from_switchpoint_cython(tree, switches_ni, ar, debug=False)[source]
ivy.chars.mk_mr.lf_mk_mr_midbranch(tree, chars, Qtype, nregime, pi=u'Equal', findmin=True)[source]
ivy.chars.mk_mr.lf_mk_mr_midbranch_mods(tree, chars, mods, pi=u'Equal', findmin=True, orderedparams=True)[source]

Create a likelihood function for testing the parameter values of user- specified models with switchpoints allowed in the middle of branches

ivy.chars.mk_mr.locs_from_switchpoint(tree, switches, locs=None)[source]

Given a tree and a list of nodes to be switchpoints, return an array of all node indices in one regime vs the other

Parameters:
  • tree (Node) – Root node of tree
  • switches (list) – List of internal nodes to act as switchpoints. Switchpoint nodes are part of their own regimes.
  • locs (np.array) – Optional pre-allocated array to store output.
Returns:

Array of indices of all nodes associated with each switchpoint, plus

all nodes “outside” the switches in the last list
Return type:

np.array
ivy.chars.mk_mr.make_modelorder_stoch(mods, name='modorder')[source]
ivy.chars.mk_mr.make_switchpoint_stoch(seg_map, name='switch')[source]
ivy.chars.mk_mr.mk_mr_midbranch(tree, chars, Qs, switchpoint, pi=u'Equal', returnPi=False, ar=None, debug=False)[source]

Calculate likelhiood of mk model with BAMM-like multiple regimes

Parameters:
  • tree (Node) – Root node of a tree. All branch lengths must be greater than 0 (except root)
  • chars (dict) –

    Dict mapping character states to tip labels. Character states should be coded 0,1,2...

    Can also be a list with tip states in preorder sequence

  • Qs (np.array) – Array of instantaneous rate matrices
  • locs (np.array) – Array of the same length as Qs containing the node indices that correspond to each Q matrix
  • p (np.array) – Optional pre-allocated p matrix
  • pi (str or np.array) –

    Option to weight the root node by given values. Either a string containing the method or an array of weights. Weights should be given in order.

    Accepted methods of weighting root:

    Equal: flat prior Equilibrium: Prior equal to stationary distribution

    of Q matrix
    Fitzjohn: Root states weighted by how well they
    explain the data at the tips.
Returns:

Log-likelihood of tree/characters given the Q matrices.
Return type:

(float)
ivy.chars.mk_mr.mk_multi_bayes(tree, chars, nregime, qidx, pi=u'Equal', seglen=0.02, stepsize=0.05)[source]

Create a Bayesian multi-mk model. User specifies which regime models to use and the Bayesian model finds the switchpoints.

Parameters:
  • tree (Node) – Root node of tree.
  • chars (dict) – Dict mapping tip labels to discrete character states. Character states must be in the form of [0,1,2...]
  • regime (int) – The number of distinct regimes to test. Set to 1 for an Mk model, set to greater than 1 for a multi-regime Mk model.
  • qidx (np.array) –

    Index specifying the model to test

    columns:
    0, 1, 2 - index axes of q 3 - index of params

    This scheme allows flexible specification of models. E.g.: Symmetric mk2:

    params = [0.2]; qidx = [[0,0,1,0],
    [0,1,0,0]]
    Asymmetric mk2:
    params = [0.2,0.6]; qidx = [[0,0,1,0],
    [0,1,0,1]]

    Note

    The qidx corresponding to the first q matrix (first column 0) is always the root regime

  • pi (str or np.array) –

    Option to weight the root node by given values. Either a string containing the method or an array of weights. Weights should be given in order.

    Accepted methods of weighting root:

    Equal: flat prior Equilibrium: Prior equal to stationary distribution

    of Q matrix
    Fitzjohn: Root states weighted by how well they
    explain the data at the tips.
  • seglen (float) – Size of segments to break tree into. The smaller this value, the more “fine-grained” the analysis will be. Optional, defaults to 2% of the root-to-tip length.
  • stepsize (float) – Maximum size of steps for switchpoints to take. Optional, defaults to 5% of root-to-tip length.
ivy.chars.mk_mr.random_tree_location(seg_map)[source]

Select random location on tree with uniform probability

Returns:
node and float, which represents the how far along the branch to
the parent node the switch occurs on
Return type:tuple
ivy.chars.mk_mr.round_step_size(step_size, seg_size)[source]

Round step_size to the nearest segment

ivy.chars.mk_mr.tree_map(tree, seglen=0.02)[source]

Make a map of the tree cut into segments of size (seglen*root_to_tip_length)

ivy.chars.recon module

ivy.chars.recon.anc_recon_cat(tree, chars, Q, p=None, pi=u'Equal', ars=None, nregime=1)[source]

Given tree, character states at tips, and transition matrix perform ancestor reconstruction for a discrete character.

Perform downpass using mk function, then perform uppass.

Return reconstructed states - including tips

Parameters:
  • tree (Node) – Root node of a tree. All branch lengths must be greater than 0 (except root)
  • chars (dict) –

    Dict mapping character states to tip labels. Character states should be coded 0,1,2...

    Can also be a list with tip states in preorder sequence

  • Q (np.array) – Instantaneous rate matrix
  • p (np.array) – 3-D array of dimensions branch_number * nchar * nchar. Optional. Pre-allocated space for efficient calculations
  • pi (str or np.array) –

    Option to weight the root node by given values. Either a string containing the method or an array of weights. Weights should be given in order. Accepted methods of weighting root:

    Equal: Flat prior Equilibrium: Prior equal to stationary distribution
    of Q matrix
    Fitzjohn: Root states weighted by how well they
    explain the data at the tips.
  • ars (dict) – Dict of pre-allocated arrays to improve speed by avoiding creating and destroying new arrays. Can be created with create_ancrecon_ars function. Optional.
  • nregime (int) – If hidden-rates model is used, the number of regimes to consider in the reconstruction. Optional, defaults to 1 (no hidden rates)
Returns:

Array of nodes in preorder sequence containing marginal

likelihoods.
Return type:

np.array
ivy.chars.recon.anc_recon_mkmr(root, chars, Q, switchpoint)[source]
ivy.chars.recon.average_nchanges(ch, pa)[source]

For use in pscore function

ivy.chars.recon.create_ancrecon_ars(tree, chars, nregime=1)[source]

Create nodelists. For use in recon function

Returns edgelist of nodes in postorder sequence, edgelist of nodes in preorder sequence, partial likelihood vector for each node’s children, childlist, and branch lengths.

ivy.chars.recon.get_childstates(node, precon, chars)[source]

For use in pscore function

ivy.chars.recon.parsimony_recon(tree, chars)[source]

Use parsimony to reconstruct the minimum number of changes needed to observe states at the tips.

Equal weights assumed

Parameters:
  • tree (Node) – Root node of tree
  • chars (dict) –

    Dict mapping character states to tip labels. Character states should be coded 0,1,2...

    Can also be a list with tip states in preorder sequence
Returns:

Array of assigned states of interior nodes
Return type:

np.array
ivy.chars.recon.pscore(tree, chars)[source]

Return the minimum number of changes needed by parsimony to observe given data on the tree

Parameters:
  • tree (Node) – Root node of the tree
  • chars (dict) – Dict mapping character states to tip labels. Character states should be coded 0,1,2...

Module contents