mvSLOUCH (version 2.7.7)

simulate_clustered_phylogeny: Simulate a phylogenetic tree with a specified number of clades.

Description

Simulate a phylogenetic tree that has a given number of clades, each with a given number of tips.

Usage

simulate_clustered_phylogeny(v_sizeclusts, joining_branchlengths = NULL, 
f_simclustphyl = "sim.bd.taxa_Yule1", joiningphyl = NULL, 
b_change_joining_branches = FALSE, ...)

Value

The resulting object is a clustered_phylo object which inherits from the phylo class and enhances it. Apart from the standard phylo fields it has two additional ones:

edges_clusters

a named list with length equalling the number of clades/clusters plus 1. The first element of the list is called joining_tree and contains the indices (row numbers of the edge matrix, indices of the edge_length vector) of the edges inside the subtree joining the clusters. Afterwords element (i+1) is named cluster_i and contains a numeric vector with the indices of the edges inside clade i.

tips_clusters

a named list with length equalling the number of clades/clusters. Each field of the list is a numeric vector containing the indices of the tips inside the clade. The names of element i of the list is cluster_i.

Arguments

v_sizeclusts

A vector with the sizes of the clades/clusters.

joining_branchlengths

Default NULL, if joiningphyl is NULL, then has to be provided. A vector of two numbers. The first element are the lengths of the branches of the cluster joining phylogeny leading to the clusters. The second element will be the lengths of the "internal" branches of the cluster joining phylogeny. If only a single number is provided, then all the branches of the joining phylogeny will have their lengths equal to this value.

f_simclustphyl

What function to use to simulate the phylogeny inside each cluster. The default value of "sim.bd.taxa_Yule1" corresponds to a pure birth tree generated by ape::rphylo(n=clade_size,birth=1,death=0), without a root branch otherwise the user should pass an object of class function and its parameters in place of the ... . The first parameter must be the number of contemporary leaves and be called n. The function has to return a valid phylo object.

joiningphyl

By what phylogeny are the clades to be joined by. Either NULL (default), a phylo object, the character string "sim.bd.taxa_Yule1" or an object of class function. If NULL, then they are joined by a caterpillar (comb/pectinate) phylogeny with the branch lengths as provided by the joining_branchlengths parameter. If it is a phylo object, then they will be joined by it. Importantly the number of tips of this phylogeny has to equal the number of clusters. If "sim.bd.taxa_Yule1", then the joining phylogeny is simulated as a pure birth tree with tips equalling the number of clusters by ape::rphylo(). If it is a function, then this is used and its parameters are passed through ... . The first parameter must be the number of contemporary leaves and be called n. The function has to return a valid phylo object.

b_change_joining_branches

Logical, if joining phylogeny (parameter joiningphyl) was provided or simulated, should its branches be changed according to what was provided in joining_branchlengths (if it was not NULL). By default FALSE and the branch lengths are not changed.

...

Parameters to be passed to user provided f_simclustphyl and joiningphyl functions. Unless one knows exactly what one is doing they should be passed by name. If there is a conflict of names, then one should pass wrapper functions around these functions where the names conflict is resolved.

Author

Krzysztof Bartoszek

References

Bartoszek K. and Vasterlund A. (2020) "Old Techniques for New Times": the RMaCzek package for producing Czekanowski's diagrams Biometrical Letters 57(2):89-118.

Examples

RNGversion(min(as.character(getRversion()),"3.6.1"))
set.seed(12345, kind = "Mersenne-Twister", normal.kind = "Inversion")
## We use a wrapper function for illustration
## a single phylo object
my_sim.bd.taxa<-function(n,...){
    ape::rphylo(n=n,...)
}
phyltree1<-simulate_clustered_phylogeny(v_sizeclusts=c(5,5,5),f_simclustphyl=my_sim.bd.taxa,
b_change_joining_branches=TRUE, joining_branchlengths=c(20,NA),joining=my_sim.bd.taxa,
birth=1,death=0)

RNGversion(min(as.character(getRversion()),"3.6.1"))
set.seed(12345, kind = "Mersenne-Twister", normal.kind = "Inversion")
## The below code should return the same tree as above
phyltree2<-simulate_clustered_phylogeny(v_sizeclusts=c(5,5,5),f_simclustphyl="sim.bd.taxa_Yule1",
b_change_joining_branches=TRUE, joining_branchlengths=c(20,NA),joining="sim.bd.taxa_Yule1")

## The resulting phylogeny is not ultrametric, if ultrametricity is required, then some procedure 
## has to be employed, e.g.
## phyltree1_u<-phytools::force.ultrametric(phyltree1, method="extend")
RNGversion(as.character(getRversion()))