if (FALSE) {
library(phyclust, quiet = TRUE)
ms()
# an ancestral tree
set.seed(1234)
(ret.ms <- ms(nsam = 3, opts = "-T -G 0.1"))
(tree.anc <- read.tree(text = ret.ms[3]))
tree.anc$tip.label <- paste("a", 1:K, sep = "")
# adjacent descendant trees to the ancestral tree
K <- 3
N <- 12
N.k <- c(3, 4, 5)
ms.dec <- NULL # a list to store trees of ms
tree.dec <- NULL # a list to store the trees in phylo class
tree.joint <- tree.anc
for(k in 1:K){
ms.dec[[k]] <- ms(N.k[k], opts = "-T -G 1.0")
tree.dec[[k]] <- read.tree(text = ms.dec[[k]][3])
tree.dec[[k]]$tip.label <- paste("d", k, ".", 1:N.k[k], sep = "")
tree.joint <- bind.tree(tree.joint, tree.dec[[k]],
where = which(tree.joint$tip.label ==
paste("a", k, sep = "")))
}
str(tree.joint)
# plot trees
par(mfrow = c(2, 3))
plot(tree.anc, main = paste("anc (", K, ")", sep = ""))
axis(1)
for(k in 1:K){
plot(tree.dec[[k]], main = paste("dec", k, " (", N.k[k], ")", sep = ""))
axis(1)
}
plot(tree.joint, main = paste("joint (", N, ")", sep = ""))
axis(1)
# use tbs option (an example from msdoc.pdf by Hudson, R.R.)
tbs.matrix <- matrix(c(3.0, 3.5, 5.0, 8.5), nrow = 2)
ret <- ms(nsam = 5, nreps = 2, opts = "-t tbs -r tbs 1000",
tbs.matrix = tbs.matrix)
print(ret)
}
Run the code above in your browser using DataLab