Hi Liam,
thank you very much for this. Yes, the function works beautifully for
what we wanted to do and the code you sent greatly facilitates
implemantation.
I was, however, wondering whether you can somehow "limit" the shifts to
happen in nodes, to avoid character reversals across a single branch.
Again, thanks a lot for the help.
Cheers,
Antigoni
On 10/11/2011 08:34 PM, Liam J. Revell wrote:
Hi Antigoni.
I have a function in the new version of "phytools" called
sim.history() that might help with this. This version is not
available yet on CRAN, but it can be downloaded from my webpage:
http://anolis.oeb.harvard.edu/~liam/R-phylogenetics/. The function
generates stochastic histories according to some rate matrix, and
stores the histories along the states at all terminal nodes.
This code might be helpful in your task. It simulates random trees,
for a given number of changes it computes the instantaneous rate
matrix that will produce that number of changes, and then it simulates
(repeatedly) under that rate until the desired number of changes is
obtained. The code fragment uses functions in "ape", "geiger", and
"phytools":
nchanges<-5 # desired number of changes
ntrees<-100 # desired number of trees
ntaxa<-100 # desired number of taxa
require(phytools); require(geiger) # required packages
mtrees<-list()
for(i in 1:ntrees){
# simulate a stochastic B-D tree using birthdeath.tree
tree<-drop.tip(birthdeath.tree(b=1,d=0,taxa.stop=(ntaxa+1)),as.character(ntaxa+1))
# compute the rate that results in (on average) nchanges changes
q<-nchanges/sum(tree$edge.length)
# simulate a stochastic history
mtree<-sim.history(tree,Q=matrix(c(-q,q,q,-q),2,2))
# loop until desired number of changes
while((length(unlist(mtree$maps))-nrow(mtree$edge))!=nchanges)
mtree<-sim.history(tree,Q=matrix(c(-q,q,q,-q),2,2))
mtrees[[i]]<-mtree
# for fun, plot the history
plotSimmap(mtrees[[i]])
}
class(mtrees)<-"multiPhylo"
# tip states for the kth tree are in mtrees[[k]]$states
Please let us know if this is useful or if you have any questions or
difficulties implementing.
All the best, Liam
--
Antigoni Kaliontzopoulou
CIBIO, Centro de Investigacao em Biodiversidade e Recursos Geneticos
Campus Agrario de Vairao, 4485-661 Vairao
PORTUGAL
Department of Ecology, Evolution, and Organismal Biology
Iowa State University, Ames,
Iowa 50011, USA
tel: +351 91 3086188
mail to: [email protected]
[email protected]
_______________________________________________
R-sig-phylo mailing list
[email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
