ah, I see. Thanks! Indeed, if I simulate trees with make.simmap or sim.history,
they work fine. But the trees I read in with read.simmap I still get the error:
OUwie(arv,randdat,model='BM1',simmap.tree=T)
Begin subplex optimization routine -- Starting value: 1
Error in which(edges[i, 5:(k + 4)] == 1) : subscript out of bounds
Which is why I thought they wouldn't work together. Looking at the trees using
print.default(), the only difference I can spot is that my imported trees don't
have a $node.states or a $states elements. Could this have something to do with
it?
Many thanks!
Code I used for simulations:
require(OUwie)
require(phytools)
require(geiger)
tree<-drop.tip(birthdeath.tree(b=1,d=0,taxa.stop=101),"101")
Q<-matrix(c(0.8,0.2,0.2,0.8),2,2, dimnames=list(c("A","B"),c("A","B")))
simtree <- sim.history(tree,Q)
simtrait<-data.frame(spp=names(simtree$states),state=simtree$states,trait=rnorm(100))
OUwie(simtree,simtrait,model='BMS',simmap.tree=T)
Q<-matrix(c(-1,0.8,0.2,0.8,-1.2,0.4,0.2,0.4,-0.6),3,3,
dimnames=list(c("A","B","C"),c("A","B","C")))
simtree <- sim.history(tree,Q)
simtrait<-data.frame(spp=names(simtree$states),state=simtree$states,trait=rnorm(100))
OUwie(simtree,simtrait,model='BMS',simmap.tree=T)
Abraços,
Rafael Maia
---
webpage: http://gozips.uakron.edu/~rm72
"A little learning is a dangerous thing; drink deep, or taste not the Pierian
spring." (A. Pope)
Graduate Student - Integrated Bioscience
University of Akron
http://gozips.uakron.edu/~shawkey/
On Mar 13, 2012, at 5:05 PM, Liam J. Revell wrote:
> Hi Rafael.
>
> I don't want to speak on behalf of the authors who are also on this list, but
> OUwie can read in modified "phylo" objects created by the phytools functions
> read.simmap & make.simmap.
>
> All the best, Liam
>
> --
> Liam J. Revell
> University of Massachusetts Boston
> web: http://faculty.umb.edu/liam.revell/
> email: liam.rev...@umb.edu
> blog: http://phytools.blogspot.com
>
> On 3/13/2012 4:35 PM, Rafael Maia wrote:
>> Speaking of OUwie, the paper just came out today:
>> Jeremy M. Beaulieu, Dwueng-Chwuan Jhwueng, Carl Boettiger & Brian C.
>> OMeara. MODELING STABILIZING SELECTION: EXPANDING THE
>> ORNSTEIN-UHLENBECK MODEL OF ADAPTIVE EVOLUTION. Evolution, Accepted
>> Article, DOI: 10.1111/j.1558-5646.2012.01619.x (link:
>> http://onlinelibrary.wiley.com/doi/10.1111/j.1558-5646.2012.01619.x/abstract
>> )
>>
>> Has anyone tried it out yet? Is there a way of converting SIMMAP trees
>> to the format to be used with the OUwie function?
>>
>> Thanks!
>>
>> Abraços,
>> Rafael Maia
>> ---
>> webpage: http://gozips.uakron.edu/~rm72
>> "A little learning is a dangerous thing; drink deep, or taste not the
>> Pierian spring." (A. Pope)
>> Graduate Student - Integrated Bioscience
>> University of Akron
>> http://gozips.uakron.edu/~shawkey/
>>
>> On Mar 12, 2012, at 2:39 PM, Liam J. Revell wrote:
>>
>>> You might want to try the OUWie package by Beaulieu & O'Meara
>>> (http://cran.r-project.org/web/packages/OUwie/index.html). I have not
>>> tried it yet, but it promises to do multi-optimum OU model fitting as
>>> well.
>>>
>>> All the best, Liam
>>>
>>> --
>>> Liam J. Revell
>>> University of Massachusetts Boston
>>> web: http://faculty.umb.edu/liam.revell/
>>> email: liam.rev...@umb.edu <mailto:liam.rev...@umb.edu>
>>> blog: http://phytools.blogspot.com
>>>
>>> On 3/12/2012 2:31 PM, Pascal Title wrote:
>>>> Hi all
>>>>
>>>> I have a dataset of 5 PC axes for a phylogeny with 52 tips and I
>>>> would like
>>>> to get the variance-covariance matrix under a global OU model of
>>>> evolution. I find that OU is strongly favored over BM, based on
>>>> fitContinuous in GEIGER.
>>>> However, I am having issues with convergence when trying to fit a hansen
>>>> model to my data. What can I do to get around this problem?
>>>> Alternatively,
>>>> is there another way to get at the multivariate VCV matrix other than
>>>> with
>>>> the OUCH package?
>>>>
>>>> The error is:
>>>> *unsuccessful convergence, code 1, see documentation for `optim'*
>>>> *Warning message:*
>>>> *In hansen(tree = ot, data = otd[varnames], regimes = otd["regimes"], :*
>>>> * unsuccessful convergence*
>>>>
>>>>
>>>> I have posted a small file that contains the following R script along
>>>> with
>>>> the data here<http://dl.dropbox.com/u/34644229/OUhansen.zip>.
>>>>
>>>> Any advice would be greatly appreciated! Thanks!
>>>>
>>>> require(ouch)
>>>> require(ape)
>>>>
>>>> tree<-read.nexus('tree.nex')
>>>> data<-read.csv('data.csv')
>>>> rownames(data)<-data$X
>>>> data$X<-NULL
>>>>
>>>> tree
>>>> head(data)
>>>>
>>>> colnames(data)->varnames #for hansen()
>>>>
>>>> ot<-ape2ouch(tree)
>>>> otd<-as(ot,"data.frame")
>>>> data$labels<-rownames(data)
>>>> otd<-merge(otd,data,by="labels",all=TRUE)
>>>> rownames(otd)<-otd$nodes
>>>> ot<-with(otd,ouchtree(nodes=nodes,ancestors=ancestors,times=times,labels=labels))
>>>> otd$regimes<-as.factor("global")
>>>> h1<-hansen(tree=ot,data=otd[varnames],regimes=otd["regimes"],sqrt.alpha=c(1,0,0,0,0,1,0,0,0,1,0,0,1,0,1),sigma=c(1,0,0,0,0,1,0,0,0,1,0,0,1,0,1))
>>>>
>>>>
>>>>
>>>
>>> _______________________________________________
>>> R-sig-phylo mailing list
>>> R-sig-phylo@r-project.org <mailto:R-sig-phylo@r-project.org>
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
>>
[[alternative HTML version deleted]]
_______________________________________________
R-sig-phylo mailing list
R-sig-phylo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
