Arguments to be passed to the constructor MixedGaussian when constructing a MGPM model with some of the default MGPM model types.

Description

Arguments to be passed to the constructor MixedGaussian when constructing a MGPM model with some of the default MGPM model types.

Usage

Args_MixedGaussian_MGPMDefaultModelTypes(omitGlobalSigmae_x = TRUE)

Arguments

Value

a list of named arguments. Currently only a named element Sigmae_x with specification depending on omitGlobalSigmae_x.

See Also

MGPMDefaultModelTypes

Arguments for the MixedGaussian constructor for scalar OU MGPM models.

Description

Arguments for the MixedGaussian constructor for scalar OU MGPM models.

Usage

Args_MixedGaussian_MGPMScalarOUType()

Value

a list.

Arguments for the MixedGaussian constructor for SURFACE OU MGPM models.

Description

Arguments for the MixedGaussian constructor for SURFACE OU MGPM models.

Usage

Args_MixedGaussian_MGPMSurfaceOUType()

Value

a list.

Latex representation of a model parameter or other found in a data.table object

Description

Latex representation of a model parameter or other found in a data.table object

Usage

FormatCellAsLatex(x)

Arguments

Value

a character string

Latex representation of a data.table with matrix and vectors in its cells

Description

Latex representation of a data.table with matrix and vectors in its cells

Usage

FormatTableAsLatex(x, argsXtable = list(), ...)

Arguments

Value

a character string representing a parseable latex text.

Examples

dt <- data.table::data.table(
   A = list(
         matrix(c(2, 0, 1.2, 3), 2, 2),
         matrix(c(2.1, 0, 1.2, 3.2, 1.3, 3.4), 3, 2)),
   b = c(2.2, 3.1))
print(FormatTableAsLatex(dt))

Class name for the scalar OU MGPM model type

Description

Class name for the scalar OU MGPM model type

Usage

MGPMScalarOUType()

Value

a character vector of one named element (ScalarOU)

Class name for the SURFACE OU MGPM model type

Description

Class name for the SURFACE OU MGPM model type

Usage

MGPMSurfaceOUType()

Value

a character vector of one named element (SURFACE)

Find the members in a list matching a member expression

Description

Find the members in a list matching a member expression

Usage

MatchListMembers(object, member, enclos = "?", q = "'", ...)

Arguments

Value

a named character vector, with names corresponding to the matched member quoted expressions (using the argument q as a quote symbol), and values corresponding to the 'enclos-ed' expressions after substituting the '?'.

See Also

PCMListMembers

Examples

model <- PCMBaseTestObjects$model_MixedGaussian_ab
MatchListMembers(model, "Sigma_x", "diag(model?[,,1L])")
MatchListMembers(model, "S.*_x", "diag(model?[,,1L])")
MatchListMembers(model, "Sigma_x", "model?[,,1L][upper.tri(model?[,,1L])]")
MatchListMembers(model, "a$Sigma_x", "model?[,,1L][upper.tri(model?[,,1L])]")

Create a multi-regime Gaussian model (MixedGaussian)

Description

Create a multi-regime Gaussian model (MixedGaussian)

Usage

MixedGaussian(
  k,
  modelTypes,
  mapping,
  className = paste0("MixedGaussian_", do.call(paste0, as.list(mapping))),
  X0 = structure(0, class = c("VectorParameter", "_Global"), description =
    "trait values at the root"),
  ...,
  Sigmae_x = structure(0, class = c("MatrixParameter", "_UpperTriangularWithDiagonal",
    "_WithNonNegativeDiagonal", "_Global"), description =
    "Upper triangular factor of the non-phylogenetic variance-covariance")
)

Arguments

Details

If X0 is not NULL it has no sense to use model-types including X0 as a parameter (e.g. use BM1 or BM3 instead of BM or BM2). Similarly if Sigmae_x is not NULL there is no meaning in using model-types including Sigmae_x as a parameter, (e.g. use OU2 or OU3 instead of OU or OU1).

Value

an object of S3 class className inheriting from MixedGaussian, GaussianPCM and PCM.

See Also

PCMTreeGetPartNames

PCMModels()

Create a template MixedGaussian object containing a regime for each model type

Description

Create a template MixedGaussian object containing a regime for each model type

Usage

MixedGaussianTemplate(mg, modelTypes = NULL)

Arguments

Value

a MixedGaussian with the same global parameter settings as for mg, the same modelTypes as modelTypes, and with a regime for each model type. The function will stop with an error if mg is not convertible to a MixedGaussian object or if there is a mismatch between the model types in mg and modelTypes.

Examples

mg <- MixedGaussianTemplate(PCMBaseTestObjects$model.ab.123.bSigmae_x)
mgTemplBMOU <- MixedGaussianTemplate(PCMBaseTestObjects$model.OU.BM)

Create a phylogenetic comparative model object

Description

This is the entry-point function for creating model objects within the PCMBase framework representing a single model-type with one or several model-regimes of this type associated with the branches of a tree. For mixed Gaussian phylogenetic models, which enable multiple model-types, use the MixedGaussian function.

Usage

PCM(
  model,
  modelTypes = class(model)[1],
  k = 1L,
  regimes = 1L,
  params = NULL,
  vecParams = NULL,
  offset = 0L,
  spec = NULL,
  ...
)

Arguments

Details

This is an S3 generic. The PCMBase package defines three methods for it:

PCM.PCM:

A default constructor for any object with a class inheriting from "PCM".

PCM.character:

A default PCM constructor from a character string specifying the type of model.

PCM.default:

A default constructor called when no other constructor is found. When called this constructor raises an error message.

Value

an object of S3 class as defined by the argument model.

See Also

MixedGaussian

Examples

# a Brownian motion model with one regime
modelBM <- PCM(model = "BM", k = 2)
# print the model
modelBM

# a BM model with two regimes
modelBM.ab <- PCM("BM", k = 2, regimes = c("a", "b"))
modelBM.ab

# print a single parameter of the model (in this case, the root value)
modelBM.ab$X0

# assign a value to this parameter (note that the brackets [] are necessary
# to preserve  the parameter attributes):
modelBM.ab$X0[] <- c(5, 2)

PCMNumTraits(modelBM)
PCMNumRegimes(modelBM)
PCMNumRegimes(modelBM.ab)

# number of numerical parameters in the model
PCMParamCount(modelBM)

# Get a vector representation of all parameters in the model
PCMParamGetShortVector(modelBM)

# Limits for the model parameters:
lowerLimit <- PCMParamLowerLimit(modelBM)
upperLimit <- PCMParamUpperLimit(modelBM)

# assign the model parameters at random: this will use uniform distribution
# with boundaries specified by PCMParamLowerLimit and PCMParamUpperLimit
# We do this in two steps:
# 1. First we generate a random vector. Note the length of the vector equals PCMParamCount(modelBM)
randomParams <- PCMParamRandomVecParams(modelBM, PCMNumTraits(modelBM), PCMNumRegimes(modelBM))
randomParams
# 2. Then we load this random vector into the model.
PCMParamLoadOrStore(modelBM, randomParams, 0, PCMNumTraits(modelBM), PCMNumRegimes(modelBM), TRUE)

print(modelBM)

PCMParamGetShortVector(modelBM)

# generate a random phylogenetic tree of 10 tips
tree <- ape::rtree(10)

#simulate the model on the tree
traitValues <- PCMSim(tree, modelBM, X0 = modelBM$X0)

# calculate the likelihood for the model parameters, given the tree and the trait values
PCMLik(traitValues, tree, modelBM)

# create a likelihood function for faster processing for this specific model.
# This function is convenient for calling in optim because it recieves and parameter
# vector instead of a model object.
likFun <- PCMCreateLikelihood(traitValues, tree, modelBM)
likFun(randomParams)

Quadratic polynomial parameters A, b, C, d, E, f for each node

Description

An S3 generic function that has to be implemented for every model class. This function is called by PCMLik.

Usage

PCMAbCdEf(
  tree,
  model,
  SE = matrix(0, PCMNumTraits(model), PCMTreeNumTips(tree)),
  metaI = PCMInfo(NULL, tree, model, verbose = verbose),
  verbose = FALSE
)

Arguments

Add a value to a list-valued attribute of a member or members matching a pattern

Description

Add a value to a list-valued attribute of a member or members matching a pattern

Usage

PCMAddToListAttribute(
  name,
  value,
  object,
  member = "",
  enclos = "?",
  spec = TRUE,
  inplace = TRUE,
  ...
)

Arguments

Value

if inplace is TRUE no value is returned. Otherwise, a modified version of object is returned.

Map a parametrization to its original form.

Description

This is an S3 generic that transforms the passed argument by applying the transformation rules for its S3 class.

This is an S3 generic. See 'PCMApplyTransformation._CholeskyFactor' for an example.

Usage

PCMApplyTransformation(o, ...)

Arguments

Details

This function returns the same object if it is not transformable.

Value

a transformed version of o.

See Also

is.Transformable

Check if the PCMBase version corresponds to a dev release

Description

Check if the PCMBase version corresponds to a dev release

Usage

PCMBaseIsADevRelease()

Value

a logical

Test objects for the PCMBase package

Description

A list containing simulated trees, trait-values and model objects for tests and examples of the PCMBase package

Usage

PCMBaseTestObjects

Format

This is a list containing the following named elements representing parameters of BM, OU and MixedGaussian models with up to three traits and up to two regimes, model objects, simulated trees with partition of their nodes in up to two parts (corresponding to the two regimes), and trait data simulated on these trees.

a.H, b.H

H matrices for OU models for regimes 'a' and 'b'.

a.Theta, b.Theta

Theta vectors for OU models for regimes 'a' and 'b'.

a.Sigma_x, b.Sigma_x

Sigma_x matrices for BM and OU models for regimes 'a' and 'b'.

a.Sigmae_x, b.Sigmae_x

Sigmae_x matrices regimes 'a' and 'b'.

a.X0, b.X0

X0 vectors for regimes 'a' and 'b'.

H

an array resulting from abind(a.H, b.H).

Theta

a matrix resulting from cbind(Theta.a, Theta.b).

Sigma_x

an array resulting from abind(a.Sigma_x, b.Sigma_x).

Sigmae_x

an array resulting from abind(a.Sigmae_x, b.Sigmae_x).

model.a.1, model.a.2, model.a.3

univariate models with a single regime for each of 3 traits.

model.a.1.Omitted_X0

same as model.a.1 but omitting X0; suitable for nesting in an MGPM model.

model.a.123, model.b.123

single-regime 3-variate models.

model.a.123.Omitted_X0

single-regime 3-variate model with omitted X0 (suitable for nesting in an MGPM.

model.a.123.Omitted_X0__bSigmae_x

same as model.a.123.Omitted_X0 but with the value of Sigmae_x copied from model.b.123.

model.a.123.Omitted_X0__Omitted_Sigmae_x

same as model.a.123 but omitting X0 and Sigmae_x.

model.b.123.Omitted_X0, model.b.123.Omitted_X0__Omitted_Sigmae_x

analogical to corresponding model.a.123...

model.ab.123

a two-regime 3-variate model.

model.ab.123.bSigmae_x

a two-regime 3-variate model having Sigmae_x from b.Sigmae_x.

model_MixedGaussian_ab

a two-regime MGPM model with a local Sigmae_x for each regime.

model_MixedGaussian_ab_globalSigmae_x

a two-regime MGPM model with a global Sigmae_x.

N

number of tips in simulated trees

tree_15_tips

a tree of 15 tips used for testing clade extraction.

tree.a

a tree with one part only (one regime)

tree.ab

a tree partitioned in two parts (two regimes)

traits.a.1

trait values simulated with model.a.1.

traits.a.123

trait values simulated with model.a.123.

traits.a.2

trait values simulated with model.a.2.

traits.a.3

trait values simulated with model.a.3.

traits.ab.123

trait values simulated with model.ab.123 on tree.ab.

tree

a tree of 5 tips used for examples.

X

3-trait data for 5 tips used together with tree for examples.

model.OU.BM

a mixed Gaussian phylogenetic model for 3 traits and an OU and BM regime used in examples.

A fixed palette of n colors

Description

A fixed palette of n colors

Usage

PCMColorPalette(
  n,
  names,
  colors = structure(hcl(h = seq(15, 375, length = n + 1), l = 65, c =
    100)[seq_len(n)], names = names)
)

Arguments

Value

A vector of character strings which can be used as color specifications by R graphics functions.

Combine all member attributes of a given name into a list

Description

Combine all member attributes of a given name into a list

Usage

PCMCombineListAttribute(object, name)

Arguments

Value

a list of attribute values

Conditional distribution of a daughter node given its parent node

Description

An S3 generic function that has to be implemented for every model class.

Usage

PCMCond(
  tree,
  model,
  r = 1,
  metaI = PCMInfo(NULL, tree, model, verbose = verbose),
  verbose = FALSE
)

Arguments

Value

an object of type specific to the type of model

Conditional distribution of a daughter node given its parent node

Description

An S3 generic function that has to be implemented for every model class.

Usage

## S3 method for class 'GaussianPCM'
PCMCond(
  tree,
  model,
  r = 1,
  metaI = PCMInfo(NULL, tree, model, verbose = verbose),
  verbose = FALSE
)

Arguments

Value

For GaussianPCM models, a named list with the following members:

Variance-covariance matrix of an OU process with optional measurement error and jump at the start

Description

Variance-covariance matrix of an OU process with optional measurement error and jump at the start

Usage

PCMCondVOU(
  H,
  Sigma,
  Sigmae = NULL,
  Sigmaj = NULL,
  xi = NULL,
  e_Ht = NULL,
  threshold.Lambda_ij = getOption("PCMBase.Threshold.Lambda_ij", 1e-08)
)

Arguments

Value

a function of one numerical argument (time) and an integer indicating the branch-index that is used to check the corresponding element in xi.

Create a likelihood function of a numerical vector parameter

Description

Create a likelihood function of a numerical vector parameter

Usage

PCMCreateLikelihood(
  X,
  tree,
  model,
  SE = matrix(0, PCMNumTraits(model), PCMTreeNumTips(tree)),
  metaI = PCMInfo(X, tree, model, SE),
  positiveValueGuard = Inf
)

Arguments

Details

It is possible to specify a function for the argument metaI. This function should have three parameters (X, tree, model) and should return a metaInfo object. (see PCMInfo).

Value

a function of a numerical vector parameter called p returning the likelihood of X given the tree and the model with parameter values specified by p.

Class names for the the default PCM and MGPM model types

Description

Utility functions returning named character vector of the model class-names for the default model types used for PCM and MixedGaussian model construction.

Usage

PCMDefaultModelTypes()

MGPMDefaultModelTypes()

Value

Both, PCMFDefaultModelTypes and MGPMDefaultModelTypes return a character string vector with named elements (A,B,C,D,E,F) defined as follows (Mitov et al. 2019a):

A.

BM (H = 0, diagonal \Sigma): BM, uncorrelated traits.

B.

BM (H = 0, symmetric \Sigma): BM, correlated traits.

C.

OU (diagonal H, diagonal \Sigma): OU, uncorrelated traits.

D.

OU (diagonal H, symmetric \Sigma): OU, correlated traits, but simple (diagonal) selection strength matrix.

E.

OU (symmetric H, symmetric \Sigma): An OU with nondiagonal symmetric H and nondiagonal symmetric \Sigma.

F.

OU (asymmetric H, symmetric \Sigma): An OU with nondiagonal asymmetric H and nondiagonal symmetric \Sigma.

The only difference between the two functions is that the model types returned by PCMFDefaultModelTypes have a global parameter X0, while the model types returned by MGPMFDefaultModelTypes have an omitted parameter X0.

References

[Mitov et al. 2019a] Mitov, V., Bartoszek, K., & Stadler, T. (2019). Automatic generation of evolutionary hypotheses using mixed Gaussian phylogenetic models. Proceedings of the National Academy of Sciences of the United States of America, 35, 201813823. http://doi.org/10.1073/pnas.1813823116

See Also

Args_MixedGaussian_MGPMDefaultModelTypes

Generate a default object of a given PCM model type or parameter type

Description

This is an S3 generic. See, e.g. 'PCMDefaultObject.MatrixParameter'.

Usage

PCMDefaultObject(spec, model, ...)

Arguments

Value

a parameter or a PCM object.

Human friendly description of a PCM

Description

Human friendly description of a PCM

Usage

PCMDescribe(model, ...)

Arguments

Details

This S3 generic function is intended to be specified for user models

Value

a character string

Describe the parameters of a PCM

Description

This is an S3 generic.

Usage

PCMDescribeParameters(model, ...)

Arguments

Value

a named list with character elements corresponding to each parameter.

Given a PCM or a parameter object, extract an analogical object for a subset of the dimensions (traits) in the original object.

Description

Given a PCM or a parameter object, extract an analogical object for a subset of the dimensions (traits) in the original object.

Usage

PCMExtractDimensions(obj, dims = seq_len(PCMNumTraits(obj)), nRepBlocks = 1L)

Arguments

Details

This is an S3 generic

Value

an object of the same class as obj with a subset of obj's dimensions multiplied nRepBlocks times.

Given a PCM or a parameter object, extract an analogical object for a subset of the regimes in the original object.

Description

Given a PCM or a parameter object, extract an analogical object for a subset of the regimes in the original object.

Usage

PCMExtractRegimes(obj, regimes = seq_len(PCMNumRegimes(obj)))

Arguments

Details

This is an S3 generic

Value

an object of the same class as obj with a subset of obj's regimes

Find the S3 method for a given PCM object or class-name and an S3 generic

Description

Find the S3 method for a given PCM object or class-name and an S3 generic

Usage

PCMFindMethod(x, method = "PCMCond")

Arguments

Value

a function object corresponding to the S3 method found or an error is raised if no such function is found for the specified object and method.

Fix a parameter in a PCM model

Description

Fix a parameter in a PCM model

Usage

PCMFixParameter(model, name)

Arguments

Value

a copy of the model with added class '_Fixed' to the class of the parameter name

Generate default model types for given PCM base-classes

Description

This function calls 'PCMListParameterizations' or 'PCMListDefaultParameterizations' and generates the corresponding 'PCMParentClasses' and 'PCMSpecify' methods in the global environment.

Usage

PCMGenerateModelTypes(
  baseTypes = list(BM = "default", OU = "default", White = "all"),
  sourceFile = NULL
)

Arguments

Value

This function has side effects only and does not return a value.

See Also

PCMListDefaultParameterizations

Generate possible parameterizations for a given type of model

Description

A parameterization of a PCM of given type, e.g. OU, is a PCM-class inheriting from this type, which imposes some restrictions or transformations of the parameters in the base-type. This function generates the S3 methods responsible for creating such parameterizations, in particular it generates the definition of the methods for the two S3 generics 'PCMParentClasses' and 'PCMSpecify' for al parameterizations specified in the 'tableParameterizations' argument.

Usage

PCMGenerateParameterizations(
  model,
  listParameterizations = PCMListParameterizations(model),
  tableParameterizations = PCMTableParameterizations(model, listParameterizations),
  env = .GlobalEnv,
  useModelClassNameForFirstRow = FALSE,
  sourceFile = NULL
)

Arguments

Value

This function does not return a value. It only has a side effect by defining S3 methods in 'env'.

Value of an attribute of an object or values for an attribute found in its members

Description

Value of an attribute of an object or values for an attribute found in its members

Usage

PCMGetAttribute(name, object, member = "", ...)

Arguments

Value

if member is an empty string, attr(object, name). Otherwise, a named list containing the value for the attribute for each member in object matched by member.

Examples

PCMGetAttribute("class", PCMBaseTestObjects$model_MixedGaussian_ab)
PCMGetAttribute(
  "dim", PCMBaseTestObjects$model_MixedGaussian_ab,
  member = "$Sigmae_x")

Get a vector of all parameters (real and discrete) describing a model on a tree including the numerical parameters of each model regime, the integer ids of the splitting nodes defining the regimes on the tree and the integer ids of the model types associated with each regime.

Description

Get a vector of all parameters (real and discrete) describing a model on a tree including the numerical parameters of each model regime, the integer ids of the splitting nodes defining the regimes on the tree and the integer ids of the model types associated with each regime.

Usage

PCMGetVecParamsRegimesAndModels(model, tree, ...)

Arguments

Details

This is an S3 generic. In the default implementation, the last entry in the returned vector is the number of numerical parameters. This is used to identify the starting positions in the vector of the first splitting node.

Value

a numeric vector concatenating the result

Meta-information about a tree and trait data associated with a PCM

Description

This function pre-processes the given tree and data in order to create meta-information used during likelihood calculation.

Usage

PCMInfo(
  X,
  tree,
  model,
  SE = matrix(0, PCMNumTraits(model), PCMTreeNumTips(tree)),
  verbose = FALSE,
  preorder = NULL,
  ...
)

Arguments

Value

a named list with the following elements:

This list is passed to PCMLik.

Likelihood of a multivariate Gaussian phylogenetic comparative model with non-interacting lineages

Description

The likelihood of a PCM represents the probability density function of observed trait values (data) at the tips of a tree given the tree and the model parameters. Seen as a function of the model parameters, the likelihood is used to fit the model to the observed trait data and the phylogenetic tree (which is typically inferred from another sort of data, such as an alignment of genetic sequences for the species at the tips of the tree). The PCMLik function provides a common interface for calculating the (log-)likelihood of different PCMs. Below we denote by N the number of tips, by M the total number of nodes in the tree including tips, internal and root node, and by k - the number of traits.

Usage

PCMLik(
  X,
  tree,
  model,
  SE = matrix(0, PCMNumTraits(model), PCMTreeNumTips(tree)),
  metaI = PCMInfo(X = X, tree = tree, model = model, SE = SE, verbose = verbose),
  log = TRUE,
  verbose = FALSE
)

Arguments

Details

For efficiency, the argument metaI can be provided explicitly, because this is not supposed to change during a model inference procedure such as likelihood maximization.

Value

a numerical value with named attributes as follows:

X0

A numerical vector of length k specifying the value at the root for which the likelihood value was calculated. If the model contains a member called X0, this vector is used; otherwise the value of X0 maximizing the likelihood for the given model parameters is calculated by maximizing the quadratic polynomial 'X0 * L_root * X0 + m_root * X0 + r_root';

error

A character string with information if a numerical or other logical error occurred during likelihood calculation.

If an error occured during likelihood calculation, the default behavior is to return NA with a non-NULL error attribute. This behavior can be changed in using global options:

"PCMBase.Value.NA"

Allows to specify a different NA value such as -Inf or -1e20 which can be used in combination with log = TRUE when using optim to maximize the log-likelihood;

"PCMBase.Errors.As.Warnings"

Setting this option to FALSE will cause any error to result in calling the stop R-base function. If not caught in a tryCatch, this will cause the inference procedure to abort at the occurence of a numerical error. By default, this option is set to TRUE, which means that getOption("PCMBase.Value.NA", as.double(NA)) is returned with an error attribute and a warning is issued.

See Also

PCMInfo PCMAbCdEf PCMLmr PCMSim PCMCond

Examples

N <- 10
tr <- PCMTree(ape::rtree(N))

model <- PCMBaseTestObjects$model_MixedGaussian_ab

PCMTreeSetPartRegimes(tr, c(`11` = 'a'), setPartition = TRUE)

set.seed(1, kind = "Mersenne-Twister", normal.kind = "Inversion")
X <- PCMSim(tr, model, X0 = rep(0, 3))

PCMLik(X, tr, model)

Calculate the likelihood of a model using the standard formula for multivariate pdf

Description

Calculate the likelihood of a model using the standard formula for multivariate pdf

Usage

PCMLikDmvNorm(
  X,
  tree,
  model,
  SE = matrix(0, PCMNumTraits(model), PCMTreeNumTips(tree)),
  metaI = PCMInfo(X, tree, model, SE, verbose = verbose),
  log = TRUE,
  verbose = FALSE
)

Arguments

Value

a numerical value with named attributes as follows:

Tracing the log-likelihood calculation of a model over each node of the tree

Description

This is an S3 generic function providing tracing information for the likelihood calculation for a given tree, data and model parameters. Useful for illustration or for debugging purpose.

Usage

PCMLikTrace(
  X,
  tree,
  model,
  SE = matrix(0, PCMNumTraits(model), PCMTreeNumTips(tree)),
  metaI = PCMInfo(X = X, tree = tree, model = model, SE = SE, verbose = verbose),
  log = TRUE,
  verbose = FALSE
)

Arguments

Value

The returned object will, in general, depend on the type of model and the algorithm used for likelihood calculation. For a G_LInv model and pruning-wise likelihood calculation, the returned object will be a data.table with columns corresponding to the node-state variables, e.g. the quadratic polynomial coefficients associated with each node in the tree.

See Also

PCMInfo PCMAbCdEf PCMLmr PCMSim PCMCond PCMParseErrorMessage

A vector of access-code strings to all members of a named list

Description

A vector of access-code strings to all members of a named list

Usage

PCMListMembers(
  l,
  recursive = TRUE,
  format = c("$", "$'", "$\"", "$`", "[['", "[[\"", "[[`")
)

Arguments

Value

a vector of character strings denoting each named member of the list.

Examples

PCMListMembers(PCMBaseTestObjects$model_MixedGaussian_ab)
PCMListMembers(PCMBaseTestObjects$model_MixedGaussian_ab, format = '$`')
PCMListMembers(PCMBaseTestObjects$tree.ab, format = '$`')

Specify the parameterizations for each parameter of a model

Description

These are S3 generics. 'PCMListParameterizations' should return all possible parametrizations for the class of 'model'. 'PCMListDefaultParameterizations' is a handy way to specify a subset of all parametrizations. 'PCMListDefaultParameterizations' should be used to avoid generating too many model parametrizations which occupy space in the R-global environment while they are not used (see PCMGenerateParameterizations). It is mandatory to implement a specification for 'PCMListParameterizations' for each newly defined class of models. 'PCMListDefaultParameterizations' has a default implementation that calls 'PCMListParameterizations' and returns the first parametrization for each parameter. Hence, implementing a method for 'PCMListDefaultParameterizations' for a newly defined model type is optional.

Usage

PCMListParameterizations(model, ...)

PCMListDefaultParameterizations(model, ...)

Arguments

Value

a named list with list elements corresponding to each parameter in model. Each list element is a list of character vectors, specifying the possible S3 class attributes for the parameter in question. For an example, type 'PCMListParameterizations.BM' to see the possible parameterizations for the BM model.

See Also

PCMGenerateParameterizations

Quadratic polynomial parameters L, m, r

Description

Quadratic polynomial parameters L, m, r

Usage

PCMLmr(
  X,
  tree,
  model,
  SE = matrix(0, PCMNumTraits(model), PCMTreeNumTips(tree)),
  metaI = PCMInfo(X = X, tree = tree, model = model, SE = SE, verbose = verbose),
  root.only = TRUE,
  verbose = FALSE
)

Arguments

Value

A list with the members A,b,C,d,E,f,L,m,r for all nodes in the tree or only for the root if root.only=TRUE.

Integer vector giving the model type index for each regime

Description

Integer vector giving the model type index for each regime

Usage

PCMMapModelTypesToRegimes(model, tree, ...)

Arguments

Details

This is a generic S3 method. The default implementation for the basic class PCM returns a vector of 1's, because it assumes that a single model type is associated with each regime. The implementation for mixed Gaussian models returns the mapping attribute of the MixedGaussian object reordered to correspond to PCMTreeGetPartNames(tree).

Value

an integer vector with elements corresponding to the elements in PCMTreeGetPartNames(tree)

Expected mean vector at each tip conditioned on a trait-value vector at the root

Description

Expected mean vector at each tip conditioned on a trait-value vector at the root

Usage

PCMMean(
  tree,
  model,
  X0 = model$X0,
  metaI = PCMInfo(NULL, tree, model, verbose = verbose),
  internal = FALSE,
  verbose = FALSE
)

Arguments

Value

If internal is FALSE (default), then a k x N matrix Mu, such that Mu[, i] equals the expected mean k-vector at tip i, conditioned on X0 and the tree. Otherwise, a k x M matrix Mu containing the mean vector for each node.

Examples

# a Brownian motion model with one regime
modelBM <- PCM(model = "BM", k = 2)
# print the model
modelBM
# assign the model parameters at random: this will use uniform distribution
# with boundaries specified by PCMParamLowerLimit and PCMParamUpperLimit
# We do this in two steps:
# 1. First we generate a random vector. Note the length of the vector equals PCMParamCount(modelBM)
randomParams <- PCMParamRandomVecParams(modelBM, PCMNumTraits(modelBM), PCMNumRegimes(modelBM))
randomParams
# 2. Then we load this random vector into the model.
PCMParamLoadOrStore(modelBM, randomParams, 0, PCMNumTraits(modelBM), PCMNumRegimes(modelBM), TRUE)

# create a random tree of 10 tips
tree <- ape::rtree(10)
PCMMean(tree, modelBM)

Calculate the mean at time t, given X0, under a PCM model

Description

Calculate the mean at time t, given X0, under a PCM model

Usage

PCMMeanAtTime(
  t,
  model,
  X0 = model$X0,
  regime = PCMRegimes(model)[1L],
  verbose = FALSE
)

Arguments

Value

A numeric vector of length k

Examples

# a Brownian motion model with one regime
modelBM <- PCM(model = "BM", k = 2)
# print the model
modelBM
# assign the model parameters at random: this will use uniform distribution
# with boundaries specified by PCMParamLowerLimit and PCMParamUpperLimit
# We do this in two steps:
# 1. First we generate a random vector. Note the length of the vector equals PCMParamCount(modelBM)
randomParams <- PCMParamRandomVecParams(modelBM, PCMNumTraits(modelBM), PCMNumRegimes(modelBM))
randomParams
# 2. Then we load this random vector into the model.
PCMParamLoadOrStore(modelBM, randomParams, 0, PCMNumTraits(modelBM), PCMNumRegimes(modelBM), TRUE)

# PCMMeanAtTime(1, modelBM)

# note that the variance at time 0 is not the 0 matrix because the model has a non-zero
# environmental deviation
PCMMeanAtTime(0, modelBM)

Get the model type(s) of a model

Description

For a regular PCM object, the model type is its S3 class. For a MixedGaussian each regime is mapped to one of several possible model types.

Usage

PCMModelTypes(obj)

Arguments

Value

a character vector

Get a list of PCM models currently implemented

Description

Get a list of PCM models currently implemented

Usage

PCMModels(pattern = NULL, parentClass = NULL, ...)

Arguments

Details

The function is using the S3 api function methods looking for all registered implementations of the function PCMSpecify.

Value

a character vector of the model classes found.

Examples

PCMModels()
PCMModels("^OU")

Number of regimes in a obj

Description

Number of regimes in a obj

Usage

PCMNumRegimes(obj)

Arguments

Value

an integer

Number of traits modeled by a PCM

Description

Number of traits modeled by a PCM

Usage

PCMNumTraits(model)

Arguments

Value

an integer

Global options for the PCMBase package

Description

Global options for the PCMBase package

Usage

PCMOptions()

Value

a named list with the currently set values of the following global options:

PCMBase.Value.NA

NA value for the likelihood; used in GaussianPCM to return this value in case of an error occurring during likelihood calculation. By default, this is set to as.double(NA).

PCMBase.Errors.As.Warnings

a logical flag indicating if errors (occuring, e.g. during likelihood calculation) should be treated as warnings and added as an attribute "error" to attach to the likelihood values. Default TRUE.

PCMBase.Raise.Lik.Errors

Should numerical and other sort of errors occurring during likelihood calculation be raised either as errors or as warnings, depending on the option PCMBase.Errors.As.Warnings. Default TRUE. This option can be useful if too frequent warnings get raised during a model fit procedure.

PCMBase.Threshold.Lambda_ij

a 0-threshold for abs(Lambda_i + Lambda_j), where Lambda_i and Lambda_j are eigenvalues of the parameter matrix H of an OU or other model. Default 1e-8. See PCMPExpxMeanExp.

PCMBase.Threshold.EV

A 0-threshold for the eigenvalues of the matrix V for a given branch. The V matrix is considered singular if it has eigenvalues smaller than PCMBase.Threshold.EV or when the ratio min(svdV)/max(svdV) is below PCMBase.Threshold.SV . Default is 1e-5. Treatment of branches with singular V matrix is defined by the option PCMBase.Skip.Singular.

PCMBase.Threshold.SV

A 0-threshold for min(svdV)/max(svdV), where svdV is the vector of singular values of the matrix V for a given branch. The V matrix is considered singular if it has eigenvalues smaller than PCMBase.Threshold.EV or when the ratio min(svdV)/max(svdV) is below PCMBase.Threshold.SV. Default is 1e-6. Treatment of branches with singular V matrix is defined by the option PCMBase.Skip.Singular.

PCMBase.Threshold.Skip.Singular

A double indicating if an internal branch of shorter length with singular matrix V should be skipped during likelihood calculation. Setting this option to a higher value, together with a TRUE value for the option PCMBase.Skip.Singular will result in tolerating some parameter values resulting in singular variance covariance matrix of the transition distribution. Default 1e-4.

PCMBase.Skip.Singular

A logical value indicating whether internal branches with singular matrix V and shorter than getOption("PCMBase.Threshold.Skip.Singular") should be skipped during likelihood calculation, adding their children L,m,r values to their parent node. Default TRUE. Note, that setting this option to FALSE may cause some models to stop working, e.g. the White model. Setting this option to FALSE will also cause errors or NA likelihood values in the case of trees with very short or 0-length branches.

PCMBase.Tolerance.Symmetric

A double specifying the tolerance in tests for symmetric matrices. Default 1e-8; see also isSymmetric.

PCMBase.Lmr.mode

An integer code specifying the parallel likelihood calculation mode.

PCMBase.ParamValue.LowerLimit

Default lower limit value for parameters, default setting is -10.0.

PCMBase.ParamValue.LowerLimit.NonNegative

Numeric (default: 0.0) indication the lower limit for parameters inheriting from class '_NonNegative'

s

PCMBase.ParamValue.LowerLimit.NonNegativeDiagonal

Default lower limit value for parameters corresponding to non-negative diagonal elements of matrices, default setting is 0.0.

PCMBase.ParamValue.UpperLimit

Default upper limit value for parameters, default setting is 10.0.

PCMBase.Transpose.Sigma_x

Should upper diagonal factors for variance-covariance rate matrices be transposed, e.g. should Sigma = t(Sigma_x) Sigma_x or, rather Sigma = Sigma_x t(Sigma_x)? Note that the two variants are not equal. The default is FALSE, meaning Sigma = Sigma_x t(Sigma_x). In this case, Sigma_x is not the actual upper Cholesky factor of Sigma, i.e. chol(Sigma) != Sigma_x. See also chol and UpperTriFactor. This option applies to parameters Sigma_x, Sigmae_x, Sigmaj_x and the measurement errors SE[,,i] for each measurement i when the argument SE is specified as a cube.

PCMBase.MaxLengthListCladePartitions

Maximum number of tree partitions returned by PCMTreeListCladePartitions. This option has the goal to interrupt the recursive search for new partitions in the case of calling PCMTreeListCladePartitions on a big tree with a small value of the maxCladeSize argument. By default this is set to Inf.

PCMBase.PCMPresentCoordinatesFun

A function with the same synopsis as PCMPresentCoordinates that can be specified in case of custom setting for the present coordinates for specific nodes of the tree. See PCMPresentCoordinates, and PCMInfo.

PCMBase.Use1DClasses

Logical indicating if 1D arithmetic operations should be used instead of multi-dimensional ones. This can speed-up computations in the case of a single trait. Currently, this feature is implemented only in the PCMBaseCpp R-package and only for some model types, such as OU and BM. Default: FALSE

PCMBase.PrintSubscript_u

Logical indicating if a subscript 'u' should be printed instead of a subscript 'x'. Used in PCMTable. Default: FALSE.

PCMBase.MaxNForGuessSigma_x

A real fraction number in the interval (0, 1) or an integer bigger than 1 controlling the number of tips to use for analytical calculation of the evolutionary rate matrix under a BM assumption. This option is used in the suggested PCMFit R-package. Default: 0.25.

PCMBase.UsePCMVarForVCV

Logical (default: FALSE) indicating if the function PCMTreeVCV should use PCMVar instead of ape's function vcv to calculate the phylogenetic variance covariance matrix under BM assumption. Note that setting this option to TRUE would slow down the function PCMTreeVCV considerably but may be more stable, particularly in the case of very big and deep trees, where previous ape's versions of the vcv function have thrown stack-overflow errors.

Examples

PCMOptions()

Create a function of time that calculates (1-exp(-lambda_ij*time))/lambda_ij for every element lambda_ij of the input matrix Lambda_ij.

Description

Create a function of time that calculates (1-exp(-lambda_ij*time))/lambda_ij for every element lambda_ij of the input matrix Lambda_ij.

Usage

PCMPExpxMeanExp(
  Lambda_ij,
  threshold.Lambda_ij = getOption("PCMBase.Threshold.Lambda_ij", 1e-08)
)

Arguments

Details

the function (1-exp(-lambda_ij*time))/lambda_ij corresponds to the product of the CDF of an exponential distribution with rate Lambda_ij multiplied by its mean value (mean waiting time).

Value

a function of time returning a matrix with entries formed from the above function or the limit, time, if |Lambda_{ij}|<=trehshold0.

Eigen-decomposition of a matrix H

Description

Eigen-decomposition of a matrix H

Usage

PCMPLambdaP_1(H)

Arguments

Details

The function fails with an error message if H is defective, that is, if its matrix of eigenvectors is computationally singular. The test for singularity is based on the rcond function.

Value

a list with elements as follows:

.

Sums of pairs of elements in a vector

Description

Sums of pairs of elements in a vector

Usage

PCMPairSums(lambda)

Arguments

Value

a squared symmetric matrix with elem_ij=lambda_i+lambda_j.

Module PCMParam

Description

Global and S3 generic functions for manipulating model parameters. The parameters in a PCM are named objects with a class attribute specifying the main type and optional properties (tags).

S3 generic functions:

PCMParamCount()

Counting the number of actual numeric parameters (used, e.g. for calculating information scores, e.g. AIC);

PCMParamLoadOrStore(), PCMParamLoadOrStore()

Storing/loading a parameter to/from a numerical vector;

PCMParamLowerLimit(),PCMParamUpperLimit()

Specifying parameter upper and lower limits;

PCMParamRandomVecParams()

Generating a random parameter vector;

For all the above properties, check-functions are defined, e.g. 'is.Local(o)', 'is.Global(o)', 'is.ScalarParameter(o)', 'is.VectorParameter', etc.

Bind named vectors or matrices into an array so that the names form the names of the last dimension.

Description

Bind named vectors or matrices into an array so that the names form the names of the last dimension.

Usage

PCMParamBindRegimeParams(...)

Arguments

Value

an array with dim attribute one longer than the number of dimensions of each argument, i.e. if there are 5 vector arguments of length 2, the returned value will be an array with dim c(2,5); if there are 5 matrix arguments of dim 2 x 2, the returned value will be an array with dim c(2,2,5).

Examples

# regimes

# in regime 'a' the three traits evolve according to three independent OU processes
a.X0 <- c(5, 2, 1)
a.H <- rbind(
  c(0, 0, 0),
  c(0, 2, 0),
  c(0, 0, 3))
a.Theta <- c(10, 6, 2)
a.Sigma_x <- rbind(
  c(1.6, 0.0, 0.0),
  c(0.0, 2.4, 0.0),
  c(0.0, 0.0, 2.0))
a.Sigmae_x <- rbind(
  c(0.0, 0.0, 0.0),
  c(0.0, 0.0, 0.0),
  c(0.0, 0.0, 0.0))
a.h_drift<-c(0, 0, 0)

# in regime 'b' there is correlation between the traits
b.X0 <- c(12, 4, 3)
b.H <- rbind(
  c(2.0, 0.1, 0.2),
  c(0.1, 0.6, 0.2),
  c(0.2, 0.2, 0.3))
b.Theta <- c(10, 6, 2)
b.Sigma_x <- rbind(
  c(1.6, 0.3, 0.3),
  c(0.0, 0.3, 0.4),
  c(0.0, 0.0, 2.0))
b.Sigmae_x <- rbind(
  c(0.2, 0.0, 0.0),
  c(0.0, 0.3, 0.0),
  c(0.0, 0.0, 0.4))
b.h_drift<-c(1, 2, 3)

H <- PCMParamBindRegimeParams(a = a.H, b = b.H)
Theta <- PCMParamBindRegimeParams(a = a.Theta, b = b.Theta)
Sigma_x <- PCMParamBindRegimeParams(a = a.Sigma_x, b = b.Sigma_x)
Sigmae_x <- PCMParamBindRegimeParams(a = a.Sigmae_x, b = b.Sigmae_x)
h_drift <- PCMParamBindRegimeParams(a = a.h_drift, b = b.h_drift)

model.a.BM_drift.123 <- PCM("BM_drift", k = 3, regimes = "a",
params = list(
  X0 = a.X0,
  h_drift = h_drift[,'a',drop=FALSE],
  Sigma_x = Sigma_x[,,'a',drop=FALSE],
  Sigmae_x = Sigmae_x[,,'a',drop=FALSE]))

# regimes 'a' and 'b', traits 1, 2 and 3
model.ab.123 <- PCM("OU", k = 3, regimes = c("a", "b"),
                    params = list(
                      X0 = a.X0,
                      H = H[,,,drop=FALSE],
                      Theta = Theta[,,drop=FALSE],
                      Sigma_x = Sigma_x[,,,drop=FALSE],
                      Sigmae_x = Sigmae_x[,,,drop=FALSE]))

Count the number of free parameters associated with a PCM or a PCM-parameter

Description

Count the number of free parameters associated with a PCM or a PCM-parameter

Usage

PCMParamCount(
  o,
  countRegimeChanges = FALSE,
  countModelTypes = FALSE,
  offset = 0L,
  k = 1L,
  R = 1L,
  parentModel = NULL
)

Arguments

Value

an integer

Get a vector of the variable numeric parameters in a model

Description

The short vector of the model parameters does not include the nodes in the tree where a regime change occurs, nor the the model types associated with each regime.

Usage

PCMParamGetShortVector(o, k = 1L, R = 1L, ...)

Arguments

Value

a numeric vector of length equal to 'PCMParamCount(o, FALSE, FALSE, 0L, k, R)'.

Load (or store) a PCM parameter from (or to) a vector of the variable parameters in a model.

Description

Load (or store) a PCM parameter from (or to) a vector of the variable parameters in a model.

Usage

PCMParamLoadOrStore(o, vecParams, offset, k, R, load, parentModel = NULL)

Arguments

Details

This S3 generic function has both, a returned value and side effects.

Value

an integer equaling the number of elements read from vecParams. In the case of type=="custom", the number of indices bigger than offset returned by the function indices(offset, k).

Locate a named parameter in the short vector representation of a model

Description

Locate a named parameter in the short vector representation of a model

Usage

PCMParamLocateInShortVector(o, accessExpr, enclos = "?")

Arguments

Value

an integer vector of length PCMParamCount(o) with NAs everywhere except at the coordinates corresponding to the parameter in question.

Examples

model <- PCM(PCMDefaultModelTypes()["D"], k = 3, regimes = c("a", "b"))
# The parameter H is a diagonal 3x3 matrix. If this matrix is considered as
# a vector the indices of its diagonal elements are 1, 5 and 9. These indices
# are indicated as the non-NA entries in the returned vector.

PCMParamLocateInShortVector(model, "$H[,,1]")
PCMParamLocateInShortVector(model, "$H[,,'a']")
PCMParamLocateInShortVector(model, "$H[,,'b']")
PCMParamLocateInShortVector(model, "$Sigma_x[,,'b']", enclos = 'diag(?)')
PCMParamLocateInShortVector(model, "$Sigma_x[,,'b']", enclos = '?[upper.tri(?)]')

The lower limit for a given model or parameter type

Description

This is an S3 generic function.

Usage

PCMParamLowerLimit(o, k, R, ...)

Arguments

Value

an object of the same S3 class as o representing a lower limit for the class.

Generate a random parameter vector for a model using uniform distribution between its lower and upper bounds.

Description

Generate a random parameter vector for a model using uniform distribution between its lower and upper bounds.

Usage

PCMParamRandomVecParams(
  o,
  k,
  R,
  n = 1L,
  argsPCMParamLowerLimit = NULL,
  argsPCMParamUpperLimit = NULL
)

Arguments

Value

a numeric matrix of dimension n x PCMParamCount(o).

See Also

PCMParamLimits PCMParamGetShortVector

Set model parameters from a named list

Description

Set model parameters from a named list

Usage

PCMParamSetByName(
  model,
  params,
  inplace = TRUE,
  replaceWholeParameters = FALSE,
  deepCopySubPCMs = FALSE,
  failIfNamesInParamsDontExist = TRUE,
  ...
)

Arguments

Value

If inplace is TRUE, the function only has a side effect of setting the parameters of the model object in the calling environment; otherwise the function returns a modified copy of the model object.

Parameter types

Description

The parameter types are divided in the following categories:

Main type

These are the "ScalarParameter", "VectorParameter" and "MatrixParameter" classes. Each model parameter must have a main type.

Scope/Omission

These are the "_Global" and "_Omitted" classes. Every parameter can be global for all regimes or local for a single regime. If not specified, local scope is assumed. In some special cases a parameter (e.g. Sigmae can be omitted from a model. This is done by adding "_Omitted" to its class attribute.

Constancy (optional)

These are the "_Fixed", "_Ones", "_Identity" and "_Zeros" classes.

Transformation (optional)

These are the "_Transformable", "_CholeskyFactor" and "_Schur" classes.

Other properties (optional)

These are the "_NonNegative", "_WithNonNegativeDiagonal", "_LowerTriangular", "_AllEqual", "_ScalarDiagonal", "_Symmetric", "_UpperTriangular", "_LowerTriangularWithDiagonal" and "_UpperTriangularWithDiagonal" classes.

Usage

is.Local(o)

is.Global(o)

is.ScalarParameter(o)

is.VectorParameter(o)

is.MatrixParameter(o)

is.WithCustomVecParams(o)

is.Fixed(o)

is.Zeros(o)

is.Ones(o)

is.Identity(o)

is.AllEqual(o)

is.NonNegative(o)

is.Diagonal(o)

is.ScalarDiagonal(o)

is.Symmetric(o)

is.UpperTriangular(o)

is.UpperTriangularWithDiagonal(o)

is.WithNonNegativeDiagonal(o)

is.LowerTriangular(o)

is.LowerTriangularWithDiagonal(o)

is.Omitted(o)

is.CholeskyFactor(o)

is.Schur(o)

is.Transformable(o)

is.Transformed(o)

is.SemiPositiveDefinite(o)

Arguments

Value

logical indicating if the object passed is from the type appearing in the function-name.

Functions

• is.Local:

• is.Global:

• is.ScalarParameter:

• is.VectorParameter:

• is.MatrixParameter:

• is.WithCustomVecParams:

• is.Fixed:

• is.Zeros:

• is.Ones:

• is.Identity:

• is.AllEqual:

• is.NonNegative:

• is.Diagonal:

• is.ScalarDiagonal:

• is.Symmetric:

• is.UpperTriangular:

• is.UpperTriangularWithDiagonal:

• is.WithNonNegativeDiagonal:

• is.LowerTriangular:

• is.LowerTriangularWithDiagonal:

• is.Omitted:

• is.CholeskyFactor:

• is.Schur:

• is.Transformable:

• is.Transformed:

• is.SemiPositiveDefinite:

The upper limit for a given model or parameter type

Description

This is an S3 generic function.

Usage

PCMParamUpperLimit(o, k, R, ...)

Arguments

Value

an object of the same S3 class as o representing an upper limit for the class.

Parent S3 classes for a model class

Description

Parent S3 classes for a model class

Usage

PCMParentClasses(model)

Arguments

Details

This S3 generic function is intended to be specified for user models. This function is called by the 'PCM.character' method to determine the parent classes for a given model class.

Value

a vector of character string denoting the names of the parent classes

Extract error information from a formatted error message.

Description

Extract error information from a formatted error message.

Usage

PCMParseErrorMessage(x)

Arguments

Value

Currently the function returns x unchanged.

A 2D Gaussian distribution density grid in the form of a ggplot object

Description

A 2D Gaussian distribution density grid in the form of a ggplot object

Usage

PCMPlotGaussianDensityGrid2D(
  mu,
  Sigma,
  xlim,
  ylim,
  xNumPoints = 100,
  yNumPoints = 100,
  ...
)

Arguments

Value

a ggplot object

A 2D sample from Gaussian distribution

Description

A 2D sample from Gaussian distribution

Usage

PCMPlotGaussianSample2D(mu, Sigma, numPoints = 1000, ...)

Arguments

Details

This function generates a random sample of numPoints 2d points using the function rmvnorm from the mvtnorm R-package. Then it produces a ggplot on the generated points.

Value

a ggplot object

Beautiful model description based on plotmath

Description

This is an S3 generic that produces a plotmath expression for its argument.

Usage

PCMPlotMath(o, roundDigits = 2, transformSigma_x = FALSE)

Arguments

Value

a character string.

Scatter plot of 2-dimensional data

Description

Scatter plot of 2-dimensional data

Usage

PCMPlotTraitData2D(
  X,
  tree,
  sizePoints = 2,
  alphaPoints = 1,
  labeledTips = NULL,
  sizeLabels = 8,
  nudgeLabels = c(0, 0),
  palette = PCMColorPalette(PCMNumRegimes(tree), PCMRegimes(tree)),
  scaleSizeWithTime = !is.ultrametric(tree),
  numTimeFacets = if (is.ultrametric(tree) || scaleSizeWithTime) 1L else 3L,
  nrowTimeFacets = 1L,
  ncolTimeFacets = numTimeFacets
)

Arguments

Value

a ggplot object

Determine which traits are present (active) on each node of the tree

Description

For every node (root, internal or tip) in tree, build a logical vector of length k with TRUE values for every present coordinate. Non-present coordinates arize from NA-values in the trait data. These can occur in two cases:

Missing measurements for some traits at some tips:

the present coordinates are FALSE for the corresponding tip and trait, but are full for all traits at all internal and root nodes.

non-existent traits for some species:

the FALSE present coordinates propagate towards the parent nodes - an internal or root node will have a present coordinate set to FALSE for a given trait, if all of its descendants have this coordinate set to FALSE.

These two cases have different effect on the likelihood calculation: missing measurements (NA) are integrated out at the parent nodes; while non-existent traits (NaN) are treated as reduced dimensionality of the vector at the parent node.

Usage

PCMPresentCoordinates(X, tree, metaI)

Arguments

Value

a k x M logical matrix. The function fails in case when all traits are NAs for some of the tips. In that case an error message is issued "PCMPresentCoordinates:: Some tips have 0 present coordinates. Consider removing these tips.".

See Also

PCMLik

Get the regimes (aka colors) of a PCM or of a PCMTree object

Description

Get the regimes (aka colors) of a PCM or of a PCMTree object

Usage

PCMRegimes(obj)

Arguments

Value

a character or an integer vector giving the regime names in the obj

Set an attribute of a named member in a PCM or other named list object

Description

Set an attribute of a named member in a PCM or other named list object

Usage

PCMSetAttribute(
  name,
  value,
  object,
  member = "",
  spec = TRUE,
  inplace = TRUE,
  ...
)

Arguments

Details

Calling this function can affect the attributes of multiple members matched by the member argument.

Value

if inplace is TRUE (default) nothing is returned. Otherwise, a modified version of object is returned.

Examples

model <- PCMBaseTestObjects$model_MixedGaussian_ab
PCMSetAttribute("class", c("MatrixParameter", "_Fixed"), model, "H")

Simulation of a phylogenetic comparative model on a tree

Description

Generate trait data on a tree according to a multivariate stochastic model with one or several regimes

Usage

PCMSim(
  tree,
  model,
  X0,
  SE = matrix(0, PCMNumTraits(model), PCMTreeNumTips(tree)),
  metaI = PCMInfo(X = NULL, tree = tree, model = model, SE = SE, verbose = verbose),
  verbose = FALSE
)

Arguments

Details

Internally, this function uses the PCMCond implementation for the given model class.

Value

numeric M x k matrix of values at all nodes of the tree, i.e. root, internal and tip, where M is the number of nodes: M=dim(tree$edge)[1]+1, with indices from 1 to N=length(tree$tip.label) corresponding to tips, N+1 corresponding to the root and bigger than N+1 corresponding to internal nodes. The function will fail in case that the length of the argument vector X0 differs from the number of traits specified in metaI$k. Error message: "PCMSim:: X0 must be of length ...".

See Also

PCMLik PCMInfo PCMCond

Examples

N <- 10
L <- 100.0
tr <- ape::stree(N)
tr$edge.length <- rep(L, N)
for(epoch in seq(1, L, by = 1.0)) {
  tr <- PCMTreeInsertSingletonsAtEpoch(tr, epoch)
}

model <- PCMBaseTestObjects$model_MixedGaussian_ab

PCMTreeSetPartRegimes(tr, c(`11` = 'a'), setPartition = TRUE)

set.seed(1, kind = "Mersenne-Twister", normal.kind = "Inversion")
X <- PCMSim(tr, model, X0 = rep(0, 3))

Parameter specification of PCM model

Description

The parameter specification of a PCM model represents a named list with an entry for each parameter of the model. Each entry in the list is a structure defining the S3 class of the parameter and its verbal description. This is an S3 generic. See 'PCMSpecify.OU' for an example method.

Usage

PCMSpecify(model, ...)

Arguments

Value

a list specifying the parameters of a PCM.

A data.table representation of a PCM object

Description

A data.table representation of a PCM object

Usage

PCMTable(
  model,
  skipGlobalRegime = FALSE,
  addTransformed = TRUE,
  removeUntransformed = TRUE
)

Arguments

Details

This is an S3 generic.

Value

an object of S3 class PCMTable

Cartesian product of possible parameterizations for the different parameters of a model

Description

This function generates a data.table in which each column corresponds to one parameter of model and each row corresponds to one combination of parameterizations for the model parameters, such that the whole table corresponds to the Cartesian product of the lists found in 'listParameterizations'. Usually, subsets of this table should be passed to 'PCMGenerateParameterizations'

Usage

PCMTableParameterizations(
  model,
  listParameterizations = PCMListParameterizations(model, ...),
  ...
)

Arguments

Value

a data.table object.

Generate a trajectory for the mean in one regime of a PCM

Description

Generate a trajectory for the mean in one regime of a PCM

Usage

PCMTrajectory(
  model,
  regime = PCMRegimes(model)[1],
  X0 = rep(0, PCMNumTraits(model)),
  W0 = matrix(0, nrow = PCMNumTraits(model), ncol = PCMNumTraits(model)),
  tX = seq(0, 100, by = 1),
  tVar = tX[seq(0, length(tX), length.out = 4)],
  dims = seq_len(PCMNumTraits(model)),
  sizeSamp = 100,
  doPlot2D = FALSE,
  plot = NULL
)

Arguments

Value

if doPlot2D is TRUE, returns a ggplot. Otherwise a named list of two elements:

dt

a data.table with columns 'regime', 't', 'X', 'V' and 'samp'. For each row corresponding to time in tVar, the column samp represents a list of sizeSamp k-vectors.

dtPlot

a data.table with the same data as in dt, but with converted columns X and samp into 2 x k columns denoted xi, i=1,...,k and xsi (i=1...k) This is suitable for plotting with ggplot.

Examples

set.seed(1, kind = "Mersenne-Twister", normal.kind = "Inversion")

# a Brownian motion model with one regime
modelOU <- PCM(model = PCMDefaultModelTypes()['F'], k = 2)

# assign the model parameters at random: this will use uniform distribution
# with boundaries specified by PCMParamLowerLimit and PCMParamUpperLimit
# We do this in two steps:
# 1. First we generate a random vector. Note the length of the vector equals
# PCMParamCount(modelBM).

randomParams <- PCMParamRandomVecParams(
  modelOU, PCMNumTraits(modelOU), PCMNumRegimes(modelOU))
# 2. Then we load this random vector into the model.
PCMParamLoadOrStore(
  modelOU,
  randomParams,
  0, PCMNumTraits(modelBM), PCMNumRegimes(modelBM), load = TRUE)

# let's plot the trajectory of the model starting from X0 = c(0,0)
PCMTrajectory(
  model = modelOU,
  X0 = c(0, 0),
  doPlot2D = TRUE)


# A faceted grid of plots for the two regimes in a mixed model:
pla <- PCMTrajectory(
  model = PCMBaseTestObjects$model_MixedGaussian_ab, regime = "a",
  X0 = c(0, 0, 0),
  doPlot2D = TRUE) +
  ggplot2::scale_y_continuous(limits = c(0, 10)) +
  ggplot2::facet_grid(.~regime)

plb <- PCMTrajectory(
  model = PCMBaseTestObjects$model_MixedGaussian_ab, regime = "b",
  X0 = c(0, 0, 0),
  doPlot2D = TRUE) +
  ggplot2::scale_y_continuous(limits = c(0, 10)) +
  ggplot2::facet_grid(.~regime) +
  ggplot2::theme(
    axis.title.y = ggplot2::element_blank(),
    axis.text.y = ggplot2::element_blank(),
    axis.ticks.y = ggplot2::element_blank())

plot(pla)
plot(plb)

Create a PCMTree object from a phylo object

Description

PCMTree is class that inherits from the class 'phylo' in the R-package 'ape'. Thus, all the functions working on a phylo object would work in the same way if they receive as argument an object of class 'PCMTree'. A PCMTree object has the following members in addition to the regular members ('tip.label', 'node.label', 'edge', 'edge.length') found in a regular phylo object:

edge.part

a character vector having as many elements as there are branches in the tree (corresponding to the rows in 'tree$edge'). Each element denotes the name of the part to which the corresponding branch belongs. A part in the tree represents a connected subset of its nodes and the branches leading to these nodes. A partition of the tree represents the splitting of the tree into a number of parts. Visually, a partition can be represented as a coloring of the tree, in which no color is assigned to more than one part. In other words, if two branches in the tree are connected by the same color, they either share a node, or all the branches on the path in the tree connecting these two branches have the same color. Formally, we define a partition of the tree as any set of nodes in the tree that includes the root. Each node in this set defines a part as the set of its descendant nodes that can be reached without traversing another partition node. We name each part by the label of its most ancestral node, that is, the node in it, which is closest to the root fo the tree. The value of edge.part for an edge in the tree is the name of the part that contains the node to which the edge is pointing.

part.regime

a named vector of size the number of parts in the tree. The names correspond to part-names whereas the values denote the ids or character names of regimes in a PCM object.

The constructor PCMTree() returns an object of call

Usage

PCMTree(tree)

Arguments

Value

an object of class PCMTree. This is a copy of the passed phylo object which is guaranteed to have node.label, edge.part and a part.regime entries set.

Examples

tree <- ape::rtree(8)

# the following four are NULLs
tree$node.label
tree$edge.part
tree$part.regime
tree$edge.regime

# In previous version regimes were assigned directly to the edges via
# tree$edge.regime. This is supported but not recommended anymore:

tree$edge.regime <- sample(
  letters[1:3], size = PCMTreeNumNodes(tree) - 1, replace = TRUE)

tree.a <- PCMTree(tree)
PCMTreeGetLabels(tree.a)
tree.a$node.label
tree.a$edge.part
tree.a$part.regime

# this is set to NULL - starting from PCMBase 1.2.9 all of the information
# for the regimes is stored in tree$edge.part and tree$part.regime.
tree.a$edge.regime

PCMTreeGetPartition(tree.a)
PCMTreeGetPartNames(tree.a)
PCMTreeGetPartRegimes(tree.a)

# let's see how the tree looks like
if(requireNamespace("ggtree"))
PCMTreePlot(tree.a) + ggtree::geom_nodelab() + ggtree::geom_tiplab()

# This is the recommended way to set a partition on the tree
PCMTreeSetPartition(tree.a, c(10, 12))

PCMTreeGetPartition(tree.a)
PCMTreeGetPartNames(tree.a)
PCMTreeGetPartRegimes(tree.a)



if(requireNamespace("ggtree"))
PCMTreePlot(tree.a) + ggtree::geom_nodelab() + ggtree::geom_tiplab()

PCMTreeGetPartsForNodes(tree.a, c(11, 15, 12))
PCMTreeGetPartsForNodes(tree.a, c("11", "15", "12"))

PCMTreeSetPartRegimes(tree.a, c(`9` = 'a', `10` = 'b', `12` = 'c'))

PCMTreeGetPartition(tree.a)
PCMTreeGetPartNames(tree.a)
PCMTreeGetPartRegimes(tree.a)

if(requireNamespace("ggtree"))
PCMTreePlot(tree.a) + ggtree::geom_nodelab() + ggtree::geom_tiplab()

Prune the tree leaving one tip for each or some of its parts

Description

Prune the tree leaving one tip for each or some of its parts

Usage

PCMTreeBackbonePartition(tree, partsToKeep = PCMTreeGetPartNames(tree))

Arguments

Value

a PCMTree object representing a pruned version of tree.

See Also

PCMTreeSetPartition

PCMTree

Examples

set.seed(1, kind = "Mersenne-Twister", normal.kind = "Inversion")
tree <- PCMTree(ape::rtree(25))

if(requireNamespace("ggtree"))
PCMTreePlot(tree) + ggtree::geom_nodelab(angle = 45) +
  ggtree::geom_tiplab(angle = 45)

backb <-  PCMTreeBackbonePartition(tree)

if(requireNamespace("ggtree"))
PCMTreePlot(backb) + ggtree::geom_nodelab(angle = 45) +
  ggtree::geom_tiplab(angle = 45)

tree2 <- PCMTreeSetPartRegimes(
  tree, c(`26` = "a", `28` = "b"), setPartition = TRUE,
  inplace = FALSE)

if(requireNamespace("ggtree"))
PCMTreePlot(tree2) + ggtree::geom_nodelab(angle = 45) +
  ggtree::geom_tiplab(angle = 45)

backb <-  PCMTreeBackbonePartition(tree2)

if(requireNamespace("ggtree"))
PCMTreePlot(backb) + ggtree::geom_nodelab(angle = 45) +
  ggtree::geom_tiplab(angle = 45)

tree3 <- PCMTreeSetPartRegimes(
  tree, c(`26` = "a", `28` = "b", `41` = "c"), setPartition = TRUE,
  inplace = FALSE)

if(requireNamespace("ggtree"))
PCMTreePlot(tree3) + ggtree::geom_nodelab(angle = 45) +
  ggtree::geom_tiplab(angle = 45)

backb <-  PCMTreeBackbonePartition(tree3)

if(requireNamespace("ggtree"))
PCMTreePlot(backb) + ggtree::geom_nodelab(angle = 45) +
  ggtree::geom_tiplab(angle = 45)

backb41 <-  PCMTreeBackbonePartition(tree3, partsToKeep = "41")

if(requireNamespace("ggtree"))
PCMTreePlot(backb41) + ggtree::geom_nodelab(angle = 45) +
  ggtree::geom_tiplab(angle = 45)

backbMoreNodes <- PCMTreeInsertSingletonsAtEpoch(
   backb, epoch = 3.7, minLength = 0.001)
PCMTreeGetPartRegimes(backbMoreNodes)

if(requireNamespace("ggtree"))
PCMTreePlot(backbMoreNodes) + ggtree::geom_nodelab(angle=45) +
  ggtree::geom_tiplab(angle=45)

backbMoreNodes <- PCMTreeInsertSingletonsAtEpoch(
   backbMoreNodes, epoch = 0.2, minLength = 0.001)
PCMTreeGetPartRegimes(backbMoreNodes)

if(requireNamespace("ggtree"))
PCMTreePlot(backbMoreNodes) + ggtree::geom_nodelab(angle=45) +
  ggtree::geom_tiplab(angle=45)

backbMoreNodes <- PCMTreeInsertSingletonsAtEpoch(
   backbMoreNodes, epoch = 1.2, minLength = 0.001)
PCMTreeGetPartRegimes(backbMoreNodes)

if(requireNamespace("ggtree"))
PCMTreePlot(backbMoreNodes) + ggtree::geom_nodelab(angle=45) +
  ggtree::geom_tiplab(angle=45)

Drop a clade from a phylogenetic tree

Description

Drop a clade from a phylogenetic tree

Usage

PCMTreeDropClade(
  tree,
  cladeRootNode,
  tableAncestors = NULL,
  X = NULL,
  returnList = !is.null(X),
  errorOnMissing = FALSE
)

Arguments

Value

If returnList is FALSE, a phylo object associated with the remaining tree after dropping the clade, otherise, a list with two named members :

tree

the phylo object associated with the remaining tree after dropping the clade

X

the submatrix of X with columns corresponding to the tips in the remaining tree

See Also

PCMTreeSpliAtNode PCMTreeExtractClade

Examples

set.seed(1, kind = "Mersenne-Twister", normal.kind = "Inversion")
tree <- PCMTree(ape::rtree(25))
PCMTreeSetPartRegimes(
  tree, c(`26`="a", `28`="b", `45`="c"), setPartition = TRUE)

if(requireNamespace("ggtree"))
PCMTreePlot(tree, palette=c(a = "red", b = "green", c = "blue")) +
  ggtree::geom_nodelab(angle = 45) + ggtree::geom_tiplab(angle = 45)

redGreenTree <- PCMTreeDropClade(tree, 45)
PCMTreeGetPartRegimes(redGreenTree)

if(requireNamespace("ggtree"))
PCMTreePlot(redGreenTree, palette=c(a = "red", b = "green", c = "blue")) +
  ggtree::geom_nodelab(angle = 45) + ggtree::geom_tiplab(angle = 45)

# we need to use the label here, because the node 29 in tree is not the same
# id in redGreenTree:
redGreenTree2 <- PCMTreeDropClade(redGreenTree, "29")

if(requireNamespace("ggtree"))
PCMTreePlot(redGreenTree2, palette=c(a = "red", b = "green", c = "blue")) +
  ggtree::geom_nodelab(angle = 45) + ggtree::geom_tiplab(angle = 45)

A data.table with time, part and regime information for the nodes in a tree

Description

A data.table with time, part and regime information for the nodes in a tree

Usage

PCMTreeDtNodes(tree)

Arguments

Value

a data.table with a row for each node in tree and columns as follows:

startNode

the starting node of each edge or NA_integer_ for the root

endNode

the end node of each edge or the root id for the root

startNodeLab

the character label for the startNode

endNodeLab

the character label for endNode

startTime

the time (distance from the root node) for the startNode or NA_real_ for the root node

endTime

the time (distance from the root node) for the endNode or NA_real_ for the root node

part

the part to which the edge belongs, i.e. the part of the endNode

regime

the regime to which the edge belongs, i.e. the regime of the part of the endNode

Examples

PCMTreeDtNodes(PCMBaseTestObjects$tree.ab)

A matrix with the begin and end time from the root for each edge in tree

Description

A matrix with the begin and end time from the root for each edge in tree

Usage

PCMTreeEdgeTimes(tree)

Arguments

Perfrorm nested extractions or drops of clades from a tree

Description

Perfrorm nested extractions or drops of clades from a tree

Usage

PCMTreeEvalNestedEDxOnTree(expr, tree)

Arguments

Value

the resulting phylo object from evaluating expr on tree.

Examples

set.seed(1, kind = "Mersenne-Twister", normal.kind = "Inversion")
tree <- PCMTree(ape::rtree(25))
PCMTreeSetPartRegimes(
  tree, c(`26`="a", `28`="b", `45`="c", `47`="d"), setPartition = TRUE)

if(requireNamespace("ggtree"))
PCMTreePlot(
  tree, palette=c(a = "red", b = "green", c = "blue", d = "magenta")) +
  ggtree::geom_nodelab(angle = 45) + ggtree::geom_tiplab(angle = 45)

bluePart <- PCMTreeEvalNestedEDxOnTree("D(E(tree,45),47)", tree)
PCMTreeGetPartRegimes(bluePart)

if(requireNamespace("ggtree"))
PCMTreePlot(
  bluePart, palette=c(a = "red", b = "green", c = "blue", d = "magenta")) +
  ggtree::geom_nodelab(angle = 45) + ggtree::geom_tiplab(angle = 45)

# Swapping the D and E calls has the same result:
bluePart2 <- PCMTreeEvalNestedEDxOnTree("E(D(tree,47),45)", tree)
PCMTreeGetPartRegimes(bluePart2)

if(requireNamespace("ggtree"))
PCMTreePlot(
  bluePart2, palette=c(a = "red", b = "green", c = "blue", d = "magenta")) +
  ggtree::geom_nodelab(angle = 45) + ggtree::geom_tiplab(angle = 45)

greenPart <- PCMTreeEvalNestedEDxOnTree("E(tree,28)", tree)

bgParts <- bluePart+greenPart


if(requireNamespace("ggtree")) {
PCMTreePlot(
  greenPart, palette=c(a = "red", b = "green", c = "blue", d = "magenta")) +
  ggtree::geom_nodelab(angle = 45) + ggtree::geom_tiplab(angle = 45)
PCMTreePlot(
  bluePart + greenPart, palette=c(a = "red", b = "green", c = "blue", d = "magenta")) +
  ggtree::geom_nodelab(angle = 45) + ggtree::geom_tiplab(angle = 45)
PCMTreePlot(
  greenPart + bluePart, palette=c(a = "red", b = "green", c = "blue", d = "magenta")) +
  ggtree::geom_nodelab(angle = 45) + ggtree::geom_tiplab(angle = 45)
PCMTreePlot(
  bgParts, palette=c(a = "red", b = "green", c = "blue", d = "magenta")) +
  ggtree::geom_nodelab(angle = 45) + ggtree::geom_tiplab(angle = 45)
}

Extract a clade from phylogenetic tree

Description

Extract a clade from phylogenetic tree

Usage

PCMTreeExtractClade(
  tree,
  cladeRootNode,
  tableAncestors = NULL,
  X = NULL,
  returnList = !is.null(X)
)

Arguments

Value

If returnList is FALSE, a phylo object associated with the clade, otherwise, a list with two named members :

tree

the phylo object associated with the clade

X

the submatrix of X with columns corresponding to the tips in the clade

See Also

PCMTreeSpliAtNode PCMTreeDropClade

Examples

set.seed(1, kind = "Mersenne-Twister", normal.kind = "Inversion")
tree <- PCMTree(ape::rtree(25))
PCMTreeSetPartRegimes(
  tree, c(`26`="a", `28`="b", `45`="c"), setPartition = TRUE)

if(requireNamespace("ggtree"))
PCMTreePlot(tree, palette=c(a = "red", b = "green", c = "blue")) +
  ggtree::geom_nodelab(angle = 45) + ggtree::geom_tiplab(angle = 45)

blueTree <- PCMTreeExtractClade(tree, 45)
PCMTreeGetPartRegimes(blueTree)

if(requireNamespace("ggtree"))
PCMTreePlot(blueTree, palette=c(a = "red", b = "green", c = "blue")) +
  ggtree::geom_nodelab(angle = 45) + ggtree::geom_tiplab(angle = 45)

# we need to use the label here, because the node 29 in tree is not the same
# id in redGreenTree:
blueTree2 <- PCMTreeDropClade(blueTree, "48")

if(requireNamespace("ggtree"))
PCMTreePlot(blueTree2, palette=c(a = "red", b = "green", c = "blue")) +
  ggtree::geom_nodelab(angle = 45) + ggtree::geom_tiplab(angle = 45)

The length of the branch leading to a node

Description

The length of the branch leading to a node

Usage

PCMTreeGetBranchLength(tree, daughterId)

Arguments

Value

a double denoting the length of the branch leading to daughterId

A vector of the daughter nodes for a given parent node id in a tree

Description

A vector of the daughter nodes for a given parent node id in a tree

Usage

PCMTreeGetDaughters(tree, parentId)

Arguments

Value

an integer vector of the direct descendants of parentId

Node labels of a tree

Description

Get the character labels of the tips, root and internal nodes in the tree (see Functions below).

Usage

PCMTreeGetLabels(tree)

PCMTreeGetRootLabel(tree)

PCMTreeGetNodeLabels(tree)

PCMTreeGetTipLabels(tree)

Arguments

Value

a character vector

Functions

• PCMTreeGetLabels: Get all labels in the order (tips,root,internal).

• PCMTreeGetRootLabel: Get the root label

• PCMTreeGetNodeLabels: Get the internal node labels

• PCMTreeGetTipLabels: Get the tip labels

The parent node id of a daughter node in a tree

Description

The parent node id of a daughter node in a tree

Usage

PCMTreeGetParent(tree, daughterId)

Arguments

Value

an integer denoting the parent of daughterId

Unique parts on a tree in the order of occurrence from the root to the tips (preorder)

Description

Unique parts on a tree in the order of occurrence from the root to the tips (preorder)

Usage

PCMTreeGetPartNames(tree)

Arguments

Value

a character vector.

Regime mapping for the parts in a tree

Description

Regime mapping for the parts in a tree

Usage

PCMTreeGetPartRegimes(tree)

Arguments

Value

a named vector with names corresponding to the part names in tree and values corresponding to regime names or ids.

Get the starting branch' nodes for each part on a tree

Description

Get the starting branch' nodes for each part on a tree

Usage

PCMTreeGetPartition(tree)

Arguments

Details

We call a starting branch the first branch from the root to the tips of a given part. A starting node is the node at which a starting branch ends.

Value

a named integer vector with elements equal to the starting nodes for each part in tree and names equal to the labels of these nodes.

See Also

PCMTreeSetPartition

Examples

set.seed(1, kind = "Mersenne-Twister", normal.kind = "Inversion")
PCMTreeGetPartition(PCMTree(ape::rtree(20)))

Get the parts of the branches leading to a set of nodes or tips

Description

Get the parts of the branches leading to a set of nodes or tips

Usage

PCMTreeGetPartsForNodes(tree, nodes = seq_len(PCMTreeNumNodes(tree)))

Arguments

Value

a character vector denoting the parts of the branches leading to the nodes, according to tree$edge.part.

Model regimes (i.e. colors) associated with the branches in a tree

Description

Model regimes (i.e. colors) associated with the branches in a tree

Usage

PCMTreeGetRegimesForEdges(tree)

Arguments

Value

a vector with entries corresponding to the rows in tree$edge denoting the regime associated with each branch in the tree. The type of the vector element is defined by the type of tree$part.regime.

Get the regimes of the branches leading to a set of nodes or tips

Description

Get the regimes of the branches leading to a set of nodes or tips

Usage

PCMTreeGetRegimesForNodes(tree, nodes = seq_len(PCMTreeNumNodes(tree)))

Arguments

Value

a character vector denoting the parts of the branches leading to the nodes, according to tree$edge.part.

Get the tips belonging to a part in a tree

Description

Get the tips belonging to a part in a tree

Usage

PCMTreeGetTipsInPart(tree, part)

Arguments

Value

an integer vector with the ids of the tips belonging to part

See Also

PCMTreeGetTipsInRegime, PCMTreeGetPartNames, PCMRegimes, PCMTreeGetPartRegimes, PCMTreeSetPartRegimes

Examples

set.seed(1, kind = "Mersenne-Twister", normal.kind = "Inversion")
tree <- ape::rtree(10)
regimes <- sample(letters[1:3], nrow(tree$edge), replace = TRUE)
PCMTreeSetRegimesForEdges(tree, regimes)

if(requireNamespace("ggtree"))
PCMTreePlot(tree) + ggtree::geom_nodelab() + ggtree::geom_tiplab()

part <- PCMTreeGetPartNames(tree)[1]
PCMTreeGetTipsInPart(tree, part)
print(part)

Get the tips belonging to a regime in a tree

Description

Get the tips belonging to a regime in a tree

Usage

PCMTreeGetTipsInRegime(tree, regime)

Arguments

Value

an integer vector with the ids of the tips belonging to regime.

See Also

PCMTreeGetTipsInPart, PCMTreeGetPartNames, PCMRegimes, PCMTreeGetPartRegimes, PCMTreeSetPartRegimes, PCMTreeGetPartition

Examples

set.seed(1, kind = "Mersenne-Twister", normal.kind = "Inversion")
tree <- ape::rtree(10)
regimes <- sample(letters[1:3], nrow(tree$edge), replace = TRUE)
PCMTreeSetRegimesForEdges(tree, regimes)

if(requireNamespace("ggtree"))
PCMTreePlot(tree) + ggtree::geom_nodelab() + ggtree::geom_tiplab()

regime <- PCMRegimes(tree)[1]
PCMTreeGetTipsInRegime(tree, regime)
print(regime)

Insert tips or singleton nodes on chosen edges

Description

Insert tips or singleton nodes on chosen edges

Usage

PCMTreeInsertSingletons(tree, nodes, positions)

PCMTreeInsertSingletonsAtEpoch(tree, epoch, minLength = 0.1)

PCMTreeInsertTipsOrSingletons(
  tree,
  nodes,
  positions = rep(0, length(nodes)),
  singleton = FALSE,
  tipBranchLengths = 0.01,
  nodeLabels = NULL,
  tipLabels = NULL
)

Arguments

Value

a modified copy of tree.

Functions

• PCMTreeInsertSingletonsAtEpoch:

• PCMTreeInsertTipsOrSingletons:

See Also

PCMTreeEdgeTimes PCMTreeLocateEpochOnBranches PCMTreeLocateMidpointsOnBranches

Examples

set.seed(1, kind = "Mersenne-Twister", normal.kind = "Inversion")
tree <- PCMTree(ape::rtree(25))
PCMTreeSetPartRegimes(
  tree, c(`26`="a", `28`="b", `45`="c", `47`="d"), setPartition = TRUE)

if(requireNamespace("ggtree"))
PCMTreePlot(
  tree,
  palette=c(a = "red", b = "green", c = "blue", d = "magenta")) +
  ggtree::geom_nodelab(angle = 45) + ggtree::geom_tiplab(angle = 45)

cbind(tree$edge, PCMTreeEdgeTimes(tree))

id47 <- PCMTreeMatchLabels(tree, "47")
length47 <- PCMTreeGetBranchLength(tree, id47)

# insert a singleton at 0.55 (root-ward) from node 47
tree <- PCMTreeInsertSingletons(tree, nodes = "47", positions = (length47/2))

if(requireNamespace("ggtree"))
PCMTreePlot(
  tree,
  palette=c(a = "red", b = "green", c = "blue", d = "magenta")) +
  ggtree::geom_nodelab(angle = 45) + ggtree::geom_tiplab(angle = 45)

# this fails, because the branch leading to node "47" is shorter now (0.55).
ggplot2::should_stop(
  tree <- PCMTreeInsertSingletons(
    tree, nodes = "47", positions= 2* length47 / 3))

# the tree is the same

if(requireNamespace("ggtree"))
PCMTreePlot(
  tree, palette=c(a = "red", b = "green", c = "blue", d = "magenta")) +
  ggtree::geom_nodelab(angle = 45) + ggtree::geom_tiplab(angle = 45)

# we can insert at a position within the edge:
tree <- PCMTreeInsertSingletons(tree, nodes = "47", positions = length47/3)

if(requireNamespace("ggtree"))
PCMTreePlot(
  tree, palette=c(a = "red", b = "green", c = "blue", d = "magenta")) +
  ggtree::geom_nodelab(angle = 45) + ggtree::geom_tiplab(angle = 45)


# Insert singletons at all branches crossing a given epoch. This will skip
# inserting singleton nodes where the resulting branches would be shorter
# than 0.1.
tree <- PCMTreeInsertSingletonsAtEpoch(tree, 2.3)

if(requireNamespace("ggtree"))
PCMTreePlot(
  tree, palette=c(a = "red", b = "green", c = "blue", d = "magenta")) +
  ggtree::geom_nodelab(angle = 45) + ggtree::geom_tiplab(angle = 45)

# Insert singletons at all branches crossing a given epoch
tree <- PCMTreeInsertSingletonsAtEpoch(tree, 2.3, minLength = 0.001)

if(requireNamespace("ggtree"))
PCMTreePlot(
  tree,
  palette=c(a = "red", b = "green", c = "blue", d = "magenta")) +
  ggtree::geom_nodelab(angle = 45) + ggtree::geom_tiplab(angle = 45)

Jumps in modeled traits associated with branches in a tree

Description

Jumps in modeled traits associated with branches in a tree

Usage

PCMTreeJumps(tree)

Arguments

Value

an integer vector of 0's and 1's with entries corresponding to the denoting if a jump took place at the beginning of a branch.

A list of all possible (including recursive) partitions of a tree

Description

A list of all possible (including recursive) partitions of a tree

Usage

PCMTreeListAllPartitions(
  tree,
  minCladeSize,
  skipNodes = character(),
  tableAncestors = NULL,
  verbose = FALSE
)

Arguments

Value

a list of integer vectors.

Examples

set.seed(1, kind = "Mersenne-Twister", normal.kind = "Inversion")
tree <- PCMTree(ape::rtree(10))

if(requireNamespace("ggtree"))
PCMTreePlot(tree) + ggtree::geom_nodelab() + ggtree::geom_tiplab()

# list of all partitions into parts of at least 4 tips
PCMTreeListAllPartitions(tree, 4)

# list of all partitions into parts of at least 3 tips
PCMTreeListAllPartitions(tree, 3)

# list all partitions into parts of at least 3 tips, excluding the partitions
# where node 16 is one of the partition nodes:
PCMTreeListAllPartitions(tree, minCladeSize = 3, skipNodes = "16")

A list of all possible clade partitions of a tree with a number of splitting nodes

Description

Each subset of nNodes distinct internal or tip nodes defines a partition of the branches of the tree into nNodes+1 blocks called parts. This function generates partitions where each part has nNodes splitting nodes and contains at least minCladeSize tips.

Usage

PCMTreeListCladePartitions(
  tree,
  nNodes,
  minCladeSize = 0,
  skipNodes = character(0),
  tableAncestors = NULL,
  countOnly = FALSE,
  verbose = FALSE
)

Arguments

Value

a list of integer nNodes-vectors. By default a full traversal of all partitions is done. It is possible to truncate the search to a limited number of partitions by setting the option PCMBase.MaxLengthListCladePartitions to a finite positive integer.

See Also

PCMOptions

A list of the descendants for each node in a tree

Description

A list of the descendants for each node in a tree

Usage

PCMTreeListDescendants(tree, tableAncestors = PCMTreeTableAncestors(tree))

Arguments

Details

This function has time and memory complexity O(M^2), where M is the number of nodes in the tree. It can take several minutes and gigabytes of memory on trees of more than 10000 tips.

Value

a list with unnamed elements in the order of nodes in the tree. Each element is an integer vector containing the descendant nodes (in increasing order) of the node identified by its index-number in the list.

A list of the path to the root from each node in a tree

Description

A list of the path to the root from each node in a tree

Usage

PCMTreeListRootPaths(tree, tableAncestors = PCMTreeTableAncestors(tree))

Arguments

Details

This function has time and memory complexity O(M^2), where M is the number of nodes in the tree. It can take several minutes and gigabytes of memory on trees of more than 10000 tips.

Value

a list with unnamed elements in the order of nodes in the tree. Each element is an integer vector containing the ancestors nodes on the path from the node (i) to the root of the tree in that order (the first element in the vector is the parent node of i and so on).

Find the crossing points of an epoch-time with each lineage of a tree

Description

Find the crossing points of an epoch-time with each lineage of a tree

Usage

PCMTreeLocateEpochOnBranches(tree, epoch)

Arguments

Value

a named list with an integer vector element "nodes" denoting the ending nodes for each branch crossing epoch and numeric vector element "positions" denoting the root-ward offset from each node in nodes.

Find the middle point of each branch longer than a threshold

Description

Find the middle point of each branch longer than a threshold

Usage

PCMTreeLocateMidpointsOnBranches(tree, threshold = 0)

Arguments

Value

a named list with an integer vector element "nodes" denoting the ending nodes for each branch crossing epoch and numeric vector element "positions" denoting the root-ward offset from each node in nodes.

Get the node numbers associated with tip- or node-labels in a tree

Description

Get the node numbers associated with tip- or node-labels in a tree

Usage

PCMTreeMatchLabels(tree, labels, stopIfNotFound = TRUE)

Arguments

Value

an integer vector giving the tip- or node- integer indices corresponding to labels. If stopIfNotFound is set to FALSE, this vector may contain NAs for the labels that were not found.

Examples

set.seed(1, kind = "Mersenne-Twister", normal.kind = "Inversion")
PCMTreeMatchLabels(PCMTree(ape::rtree(20)), c("t1", "t15", "21", "39"))
PCMTreeMatchLabels(PCMTree(ape::rtree(20)), c("t1", "45"), stopIfNotFound = FALSE)

Which couples from a given set of nodes in a tree belong to the same part?

Description

Which couples from a given set of nodes in a tree belong to the same part?

Which couples from a given set of nodes in a tree belong to the same regime?

Usage

PCMTreeMatrixNodesInSamePart(
  tree,
  nodes = seq_len(PCMTreeNumNodes(tree)),
  upperTriangle = TRUE,
  returnVector = TRUE
)

PCMTreeMatrixNodesInSameRegime(
  tree,
  nodes = seq_len(PCMTreeNumNodes(tree)),
  upperTriangle = TRUE,
  returnVector = TRUE
)

Arguments

Value

a L x L logical matrix with TRUE on the diagonal and for each couple of tips that belong to the same part or regime. If returnVector is TRUE (default) only a vector of the non-NA entries will be returned.

a L x L logical matrix with TRUE on the diagonal and for each couple of tips that belong to the same part or regime. If returnVector is TRUE (default) only a vector of the non-NA entries will be returned.

Examples

set.seed(1, kind = "Mersenne-Twister", normal.kind = "Inversion")
tree <- PCMTree(ape::rtree(8))
PCMTreeMatrixNodesInSamePart(tree, returnVector = FALSE)

PCMTreeSetPartition(tree, c(10, 12))
PCMTreeMatrixNodesInSamePart(tree, returnVector = FALSE)

PCMTreeMatrixNodesInSamePart(tree)
PCMTreeMatrixNodesInSamePart(tree, seq_len(PCMTreeNumTips(tree)))
PCMTreeMatrixNodesInSamePart(
  tree, seq_len(PCMTreeNumTips(tree)), returnVector = FALSE)

set.seed(1, kind = "Mersenne-Twister", normal.kind = "Inversion")
tree <- PCMTree(ape::rtree(8))
PCMTreeMatrixNodesInSamePart(tree, returnVector = FALSE)

PCMTreeSetPartition(tree, c(10, 12))
PCMTreeMatrixNodesInSamePart(tree, returnVector = FALSE)

PCMTreeMatrixNodesInSamePart(tree)
PCMTreeMatrixNodesInSamePart(tree, seq_len(PCMTreeNumTips(tree)))
PCMTreeMatrixNodesInSamePart(
  tree, seq_len(PCMTreeNumTips(tree)), returnVector = FALSE)

Find the nearest node to a given time from the root (epoch) on each lineage crossing this epoch

Description

Find the nearest node to a given time from the root (epoch) on each lineage crossing this epoch

Usage

PCMTreeNearestNodesToEpoch(tree, epoch)

Arguments

Value

an integer vector

Calculate the time from the root to each node of the tree

Description

Calculate the time from the root to each node of the tree

Usage

PCMTreeNodeTimes(tree, tipsOnly = FALSE)

Arguments

Value

A vector of size the number of nodes in the tree (tips, root, internal) containing the time from the root to the corresponding node in the tree.

Number of all nodes in a tree

Description

Number of all nodes in a tree

Usage

PCMTreeNumNodes(tree)

Arguments

Details

Wrapper for nrow(tree$edge) + 1

Value

the number of nodes in tree including root, internal and tips.

Number of unique parts on a tree

Description

Number of unique parts on a tree

Usage

PCMTreeNumParts(tree)

Arguments

Value

the number of different parts encountered on the tree branches

Wrapper for length(tree$tip.label)

Description

Wrapper for length(tree$tip.label)

Usage

PCMTreeNumTips(tree)

Arguments

Value

the number of tips in tree

Plot a tree with parts and regimes assigned to these parts

Description

Plot a tree with parts and regimes assigned to these parts

Usage

PCMTreePlot(
  tree,
  palette = PCMColorPalette(PCMNumRegimes(tree), PCMRegimes(tree)),
  ...
)

Arguments

Note

This function requires that the ggtree package is installed. At the time of releasing this version the ggtree package is not available on CRAN. Check the ggtree homepage for instruction on how to install this package: .

Post-order tree traversal

Description

Post-order tree traversal

Usage

PCMTreePostorder(tree)

Arguments

Value

a vector of indices of edges in tree$edge in post-order.

Pre-order tree traversal

Description

Pre-order tree traversal

Usage

PCMTreePreorder(tree)

Arguments

Value

a vector of indices of edges in tree$edge in pre-order.

Set tip and internal node labels in a tree

Description

Set tip and internal node labels in a tree

Usage

PCMTreeSetLabels(
  tree,
  labels = as.character(1:PCMTreeNumNodes(tree)),
  inplace = TRUE
)

Arguments

Value

if inplace is FALSE, a copy of tree with set or modified tree$tip.label and tree$node.label. If the original tree has a member edge.part, the returned tree has tree$edge.part and tree$part.regime updated. If inplace is TRUE (the default), nothing is returned and the above changes are made directly on the input tree.

See Also

PCMTree

Examples

tree <- ape::rtree(5)
tree$tip.label
# the following three are NULLs
tree$node.label
tree$edge.part
tree$part.regime


tree <- PCMTree(tree)
PCMTreeSetPartition(tree, c(6, 8))
tree$tip.label
tree$node.label
tree$edge.part
tree$part.regime

PCMTreeSetLabels(
  tree, labels = paste0(c(rep("t", 5), rep("n", 4)), PCMTreeGetLabels(tree)))
PCMTreeGetLabels(tree)
tree$tip.label
tree$node.label
tree$edge.part
tree$part.regime

Set regimes for the parts in a tree

Description

Set regimes for the parts in a tree

Usage

PCMTreeSetPartRegimes(tree, part.regime, setPartition = FALSE, inplace = TRUE)

Arguments

Value

If inplace is TRUE nothing, otherwise a copy of the tree with set edge.part and part.regime members.

See Also

PCMTree

Examples

tree <- PCMTree(ape::rtree(25))
PCMTreeGetPartition(tree)
PCMTreeGetPartRegimes(tree)
PCMTreeGetPartNames(tree)

PCMTreeSetPartRegimes(tree, c(`26` = 2))
PCMTreeGetPartition(tree)
PCMTreeGetPartRegimes(tree)
PCMTreeGetPartNames(tree)

PCMTreeSetPartRegimes(tree, c(`26` = "global-regime"))
PCMTreeGetPartition(tree)
PCMTreeGetPartRegimes(tree)
PCMTreeGetPartNames(tree)

# This should fail because no partition with nodes 26, 28 and 41 has been
# done.
ggplot2::should_stop(
  PCMTreeSetPartRegimes(tree, c(`26` = "a", `28` = "b", `41` = "c")))
# This should succeed and change the partition as well as regime assignment
PCMTreeSetPartRegimes(
  tree, c(`26` = "a", `28` = "b", `41` = "c"), setPartition = TRUE)
PCMTreeGetPartition(tree)
PCMTreeGetPartRegimes(tree)
PCMTreeGetPartNames(tree)



set.seed(1, kind = "Mersenne-Twister", normal.kind = "Inversion")
# number of tips
N <- 40

# tree with one regime
tree.a <- ape::rtree(N)

tree.a <- PCMTree(tree.a)

PCMTreeSetPartRegimes(
  tree.a,
  part.regime = structure("a", names = as.character(N+1L)),
  setPartition = TRUE,
  inplace = TRUE)



if(requireNamespace("ggtree"))
PCMTreePlot(tree.a) + ggtree::geom_nodelab() + ggtree::geom_tiplab()

tree.ab <- tree.a
PCMTreeSetPartRegimes(
  tree.ab,
  part.regime = structure(c("a", "b"), names = as.character(c(N+1L, N+31L))),
  setPartition = TRUE,
  inplace = TRUE)

if(requireNamespace("ggtree"))
PCMTreePlot(tree.ab) + ggtree::geom_nodelab() + ggtree::geom_tiplab()

Set a partition of a tree by specifying the partition nodes

Description

Set a partition of a tree by specifying the partition nodes

Usage

PCMTreeSetPartition(tree, nodes = c(PCMTreeNumTips(tree) + 1L), inplace = TRUE)

Arguments

Value

If inplace is TRUE nothing, otherwise a copy of the tree with set edge.part member.

See Also

PCMTreeGetPartition

PCMTree

Examples

set.seed(1, kind = "Mersenne-Twister", normal.kind = "Inversion")
tree <- PCMTree(ape::rtree(8))
PCMTreeSetLabels(tree, paste0("x", PCMTreeGetLabels(tree)))
PCMTreeGetPartition(tree)
PCMTreeGetPartNames(tree)
PCMTreeGetPartRegimes(tree)

if(requireNamespace("ggtree"))
PCMTreePlot(tree) + ggtree::geom_nodelab() + ggtree::geom_tiplab()

tree <- PCMTreeSetPartition(tree, c(12, 14), inplace = FALSE)
PCMTreeGetPartition(tree)
PCMTreeGetPartNames(tree)
PCMTreeGetPartRegimes(tree)

if(requireNamespace("ggtree"))
PCMTreePlot(tree) + ggtree::geom_nodelab() + ggtree::geom_tiplab()

# reset the partition to a default one, where there is only one part:
PCMTreeSetPartition(tree)

PCMTreeGetPartition(tree)
PCMTreeGetPartNames(tree)
PCMTreeGetPartRegimes(tree)
if(requireNamespace("ggtree"))
PCMTreePlot(tree) + ggtree::geom_nodelab() + ggtree::geom_tiplab()


# reset the labels to the default labels which are character representations
# of the node ids
PCMTreeSetLabels(tree)
PCMTreeGetPartition(tree)
PCMTreeGetPartNames(tree)
PCMTreeGetPartRegimes(tree)

Set the regime for each individual edge in a tree explicitly

Description

Set the regime for each individual edge in a tree explicitly

Usage

PCMTreeSetRegimesForEdges(tree, regimes, inplace = TRUE)

Arguments

Value

if inplace is TRUE, nothing, otherwise a modified copy of tree.

Note

Calling this function overwrites the current partitioning of the tree.

Examples

set.seed(1, kind = "Mersenne-Twister", normal.kind = "Inversion")
tree <- ape::rtree(10)
regimes <- sample(letters[1:3], nrow(tree$edge), replace = TRUE)
PCMTreeSetRegimesForEdges(tree, regimes)

if(requireNamespace("ggtree"))
PCMTreePlot(tree)

Slit a tree at a given internal node into a clade rooted at this node and the remaining tree after dropping this clade

Description

Slit a tree at a given internal node into a clade rooted at this node and the remaining tree after dropping this clade

Usage

PCMTreeSplitAtNode(
  tree,
  node,
  tableAncestors = PCMTreeTableAncestors(tree),
  X = NULL
)

Arguments

Details

In the current implementation, the edge.jump and edge.part members of the tree will be discarded and not present in the clade.

Value

A list containing two named phylo objects:

clade

The subtree (clade) starting at node.

Xclade

The portion of X attributable to the tips in clade; NULL if X is NULL.

rest

The tree resulting after dropping all tips in the clade.

Xrest

The portion of X attributable to the tips in rest; NULL if X is NULL.

Examples

set.seed(1, kind = "Mersenne-Twister", normal.kind = "Inversion")
tree <- PCMTree(ape::rtree(25))

if(requireNamespace("ggtree"))
PCMTreePlot(tree) + ggtree::geom_nodelab(angle = 45) +
  ggtree::geom_tiplab(angle = 45)

spl <- PCMTreeSplitAtNode(tree, 28)

if(requireNamespace("ggtree"))
PCMTreePlot(PCMTree(spl$clade)) + ggtree::geom_nodelab(angle = 45) +
  ggtree::geom_tiplab(angle = 45)

if(requireNamespace("ggtree"))
PCMTreePlot(PCMTree(spl$rest)) + ggtree::geom_nodelab(angle = 45) +
  ggtree::geom_tiplab(angle = 45)

A matrix (table) of ancestors/descendants for each node in a tree

Description

A matrix (table) of ancestors/descendants for each node in a tree

Usage

PCMTreeTableAncestors(tree, preorder = PCMTreePreorder(tree))

Arguments

Details

This function has time and memory complexity O(M^2), where M is the number of nodes in the tree. It can take several minutes and gigabytes of memory on trees of more than 10000 tips.

Value

an integer square matrix of size M x M where M is the number of nodes in the tree. Element j on row i is 0 if j is not an ancestor of i or a positive integer equal to the position of j on the path from the root to i if j is an ancestor of i.

A character representation of a phylo object.

Description

A character representation of a phylo object.

Usage

PCMTreeToString(tree, includeLengths = FALSE, includePartition = FALSE)

Arguments

Value

a character string.

Phylogenetic Variance-covariance matrix

Description

This is a simplified wrapper for ape's vcv function. Setting the runtime option PCMBase.UsePCMVarForVCV to TRUE will switch to a computation of the matrix using the function PCMVar.

Usage

PCMTreeVCV(tree)

Arguments

Value

a N x N matrix. Assuming a BM model of evolution, this is a matrix in which element (i,j) is equal to the shared root-distance of the nodes i and j.

See Also

vcv PCMVar PCMOptions

Unfix a parameter in a PCM model

Description

Unfix a parameter in a PCM model

Usage

PCMUnfixParameter(model, name)

Arguments

Value

a copy of the model with removed class '_Fixed' from the class of the parameter name

Expected variance-covariance matrix for each couple of tips (i,j)

Description

Expected variance-covariance matrix for each couple of tips (i,j)

Usage

PCMVar(
  tree,
  model,
  W0 = matrix(0, PCMNumTraits(model), PCMNumTraits(model)),
  SE = matrix(0, PCMNumTraits(model), PCMTreeNumTips(tree)),
  metaI = PCMInfo(NULL, tree, model, verbose = verbose),
  internal = FALSE,
  diagOnly = FALSE,
  verbose = FALSE
)

Arguments

Value

If internal is FALSE, a (k x N) x (k x N) matrix W, such that k x k block W[((i-1)*k)+(1:k), ((j-1)*k)+(1:k)] equals the expected covariance matrix between tips i and j. Otherwise, a list with an element 'W' as described above and a k x M matrix element 'Wii' containing the per-node variance covariance matrix for each node: The k x k block Wii[, (i-1)*k + (1:k)] represents the variance covariance matrix for node i.

Examples

# a Brownian motion model with one regime
modelBM <- PCM(model = "BM", k = 2)
# print the model
modelBM
# assign the model parameters at random: this will use uniform distribution
# with boundaries specified by PCMParamLowerLimit and PCMParamUpperLimit
# We do this in two steps:
# 1. First we generate a random vector. Note the length of the vector equals PCMParamCount(modelBM)
randomParams <- PCMParamRandomVecParams(modelBM, PCMNumTraits(modelBM), PCMNumRegimes(modelBM))
randomParams
# 2. Then we load this random vector into the model.
PCMParamLoadOrStore(modelBM, randomParams, 0, PCMNumTraits(modelBM), PCMNumRegimes(modelBM), TRUE)

# create a random tree of 10 tips
tree <- ape::rtree(10)
covMat <- PCMVar(tree, modelBM)

Calculate the variance covariance k x k matrix at time t, under a PCM model

Description

Calculate the variance covariance k x k matrix at time t, under a PCM model

Usage

PCMVarAtTime(
  t,
  model,
  W0 = matrix(0, PCMNumTraits(model), PCMNumTraits(model)),
  SE = matrix(0, PCMNumTraits(model), PCMNumTraits(model)),
  regime = PCMRegimes(model)[1L],
  verbose = FALSE
)

Arguments

Value

A numeric k x k matrix

Examples

# a Brownian motion model with one regime
modelBM <- PCM(model = "BM", k = 2)
# print the model
modelBM
# assign the model parameters at random: this will use uniform distribution
# with boundaries specified by PCMParamLowerLimit and PCMParamUpperLimit
# We do this in two steps:
# 1. First we generate a random vector. Note the length of the vector equals PCMParamCount(modelBM)
randomParams <- PCMParamRandomVecParams(modelBM, PCMNumTraits(modelBM), PCMNumRegimes(modelBM))
randomParams
# 2. Then we load this random vector into the model.
PCMParamLoadOrStore(modelBM, randomParams, 0, PCMNumTraits(modelBM), PCMNumRegimes(modelBM), TRUE)

# PCMVarAtTime(1, modelBM)

# note that the variance at time 0 is not the 0 matrix because the model has a non-zero
# environmental deviation
PCMVarAtTime(0, modelBM)

Check if all packages listed in Suggests are available

Description

Check if all packages listed in Suggests are available

Usage

RequireSuggestedPackages()

Value

logical TRUE if suggested packages are installed and can be loaded; FALSE otherwise

True positive rate of a set of binary predictions against their trues

Description

Let the set of predictions be described by a logical vector 'pred', and let the corresponding trues by described in a logical vector 'true' of the same length. Then, the true positive rate is given by the expression: sum(pred & true)/sum(true). The false positive rate is given by the expression: sum(pred & !true)/sum(!true). If these expressions do not give a finite number, NA_real_ is returned.

Usage

TruePositiveRate(pred, true)

FalsePositiveRate(pred, true)

Arguments

Value

a double between 0 and 1 or NA_real_ if the result is not a finite number.

Examples

TruePositiveRate(c(1,0,1,1), c(1,1,0,1))
TruePositiveRate(c(0,0,0,0), c(1,1,0,1))
TruePositiveRate(c(1,1,1,1), c(1,1,0,1))
FalsePositiveRate(c(1,0,1,1), c(1,1,0,1))
FalsePositiveRate(c(0,0,0,0), c(1,1,0,1))
FalsePositiveRate(c(1,1,1,1), c(1,1,0,1))
TruePositiveRate(c(1,0,1,1), c(0,0,0,0))
FalsePositiveRate(c(1,0,1,1), c(1,1,1,1))

Upper triangular factor of a symmetric positive definite matrix

Description

This function is an analog to the Cholesky decomposition.

Usage

UpperTriFactor(Sigma)

Arguments

Value

an upper triangular matrix Sigma_x, such that Sigma = Sigma_x %*% t(Sigma_x)

See Also

chol;

the option PCMBase.Transpose.Sigma_x in PCMOptions.

Examples

# S is a symmetric positive definite matrix
M<-matrix(rexp(9),3,3); S <- M %*% t(M)

# This should return a zero matrix:
UpperTriFactor(S) %*% t(UpperTriFactor(S)) - S

# This should return a zero matrix too:
t(chol(S)) %*% chol(S) - S

# Unless S is diagonal, in the general case, this will return a
# non-zero matrix:
chol(S) %*% t(chol(S)) - S

White Gaussian PCM ignoring phylogenetic history

Description

White model ignoring phylogenetic history, treating trait values as independent samples from a k-variate Gaussian.

Details

Calculating likelihoods for this model does not work if the global option PCMBase.Singular.Skip is set to FALSE.

Convert a GaussianPCM model object to a MixedGaussian model object

Description

Convert a GaussianPCM model object to a MixedGaussian model object

Usage

as.MixedGaussian(o, modelTypes = NULL)

Arguments

Value

a MixedGaussian object.

Examples

mg <- as.MixedGaussian(PCMBaseTestObjects$model.ab.123.bSigmae_x)
stopifnot(
  PCMLik(
    X = PCMBaseTestObjects$traits.ab.123,
    PCMBaseTestObjects$tree.ab,
    PCMBaseTestObjects$model.ab.123.bSigmae_x) ==
  PCMLik(
    X = PCMBaseTestObjects$traits.ab.123,
    PCMBaseTestObjects$tree.ab,
    mg))

Data for Fig3 in the TPB manuscript

Description

A list containing simulated tree, models and data used in Fig. 3

Usage

dataFig3

Format

This is a list containing the following named elements representing simulation parameters, a simulated tree and PCM objects, used in Fig. 3. For details on all these objects, read the file data-raw/Fig3.Rmd.

Check if an object is a 'GaussianPCM'

Description

Check if an object is a 'GaussianPCM'

Usage

is.GaussianPCM(x)

Arguments

Value

TRUE if x inherits from the S3 class 'GaussianPCM', FALSE otherwise.

Check if an object is a 'MixedGaussian' PCM

Description

Check if an object is a 'MixedGaussian' PCM

Usage

is.MixedGaussian(x)

Arguments

Value

TRUE if x inherits from the S3 class 'MixedGaussian', FALSE otherwise.

Check if an object is a PCM.

Description

Check if an object is a PCM.

Usage

is.PCM(x)

Arguments

Value

TRUE if 'x' inherits from the S3 class "PCM".

Check that a tree is a PCMTree

Description

Check that a tree is a PCMTree

Usage

is.PCMTree(tree)

Arguments

Value

a logical TRUE if 'inherits(tree, "PCMTree")' is TRUE.