set.seed(42)
## A simple example:
## A random tree with 10 tips
tree <- rcoal(10)
## Setting up the parameters
my_rates = c(rgamma, rate = 10, shape = 5)
## A random Mk matrix (10*50)
matrix_simple <- sim.morpho(tree, characters = 50, model = "ER", rates = my_rates,
invariant = FALSE)
## Run a basic ancestral states estimations
ancestral_states <- multi.ace(matrix_simple, tree)
ancestral_states[1:5, 1:5]
## A more complex example
## Create a multiple list of 5 trees
multiple_trees <- rmtree(5, 10)
## Modify the matrix to contain missing and special data
matrix_complex <- matrix_simple
matrix_complex[sample(1:length(matrix_complex), 50)] <- "-"
matrix_complex[sample(1:length(matrix_complex), 50)] <- "0%2"
matrix_complex[sample(1:length(matrix_complex), 50)] <- "?"
matrix_complex[1:5,1:5]
## Set a list of extra special tokens
my_spec_tokens <- c("weirdtoken" = "%")
## Set some special behaviours for the "weirdtoken" and for "-" and "?"
my_spec_behaviours <- list()
## Inapplicable tokens "-" are ignored
my_spec_behaviours$inapplicable <- function(x,y) return(NA)
## Missing tokens "?" are considered as all states
my_spec_behaviours$missing <- function(x,y) return(y)
## Weird tokens are considered as state 0 and 3
my_spec_behaviours$weirdtoken <- function(x,y) return(c(1,2))
## Create a random branch length modifier to apply to each tree
branch_lengths <- rnorm(18)^2
## Setting a list of model ("ER" for the 25 first characters and then "SYM")
my_models <- c(rep("ER", 25), rep("SYM", 25))
## Run the ancestral states on all the tree with multiple options
ancestral_states <- multi.ace(matrix_complex, multiple_trees,
verbose = TRUE,
models = my_models,
threshold = 0.95,
special.tokens = my_spec_tokens,
special.behaviours = my_spec_behaviours,
brlen.multiplier = branch_lengths,
output = "combined.matrix")
## The results for the the two first characters for the first tree
ancestral_states[[1]][, 1:2]
if (FALSE) {
## The same example but running in parallel
ancestral_states <- multi.ace(matrix_complex, multiple_trees,
verbose = TRUE,
models = my_models,
threshold = 0.95,
special.tokens = my_spec_tokens,
special.behaviours = my_spec_behaviours,
brlen.multiplier = branch_lengths,
output = "combined.matrix",
parallel = TRUE)
}
Run the code above in your browser using DataLab