Hi Agus.
It looks like in your example the Hessian matrix of partial second
derivatives (the curvature) at the optimum is numerically singular. I'm
not sure what this means - it could be that the surface is very flat, or
alternatively that it is highly correlated.
Is the error fatal, or do you get a result? The Hessian matrix is only
used to estimate standard errors of the fitted parameter values, so if
the option to compute the variances from the Hessian could be turned off
(not possible in the present version), presumably the analysis would run
through.
Another option is to use brownieREML in phytools or, as Brian pointed
out yesterday, the OUwie package. brownieREML is a "lighter" REML
version of brownie.lite, which does not compute the Hessian or the
standard errors on parameter estimates. OUwie is a separate package that
uses phytools mapped trees to do the analyses of brownie.lite in one of
its functions.
All the best, Liam
Liam J. Revell, Assistant Professor of Biology
University of Massachusetts Boston
web: http://faculty.umb.edu/liam.revell/
email: [email protected]
blog: http://phytools.blogspot.com
On 9/18/2012 4:08 PM, Agus Camacho wrote:
Dear all, I tried to test something similar to Jason, but relating the
change in rates to the acquisition of a trait, instead to a specific date.
For that, I used a phylogeny with ten taxa, without outgroup.
I used exactly the same script gently proposed by Liam but got the
following error message:
"Error in solve.default(model1$hessian) :Lapack routine dgesv: system is
exactly singular"
I've been looking through the internet but was not able to identify what
causes this error. Any hint about that?
BTW, anybody knows how to directly answer a message found in the digest of
the mail list, without leaving the thread?
Below, I pasted the complete thread dealing with this topic.
Cheers,
Agus
Message: 9
Date: Tue, 18 Sep 2012 00:13:00 -0400
From: "Brian O'Meara" <[email protected]>
To: Jason S <[email protected]>
Cc: "[email protected]" <[email protected]>
Subject: Re: [R-sig-phylo] variation in rates over time
Message-ID:
<cakywhkq-jjroqzdvfutphi4xopxam+aaacu64nldb3inxoz...@mail.gmail.com>
Content-Type: text/plain
I agree with the suggestions so far. I just wanted to point out a few more
alternatives:
You could use the geiger package to estimate the best scaling for the
tworatetree transformation to do this (should be equivalent to the earlier
solutions, though it would require running optimization).
You could also use the OUwie package with a tree that has been given simmap
mappings using phytools. The advantage of this is that you could evaluate
Brownian models but you could also look at various Ornstein-Uhlenbeck
models (though note that while there is information about different
Brownian rates before and after a time slice, info about different alphas
(strength of pull parameter) and thetas (the attractor, aka "optimal
value") is rapidly (but not immediately) lost under an OU process).
For completeness, especially for citations, note that O'Meara et al. (2006)
and Thomas et al. (2006) independently arrived at essentially the same
method, so it is worth reading both papers.
Best,
Brian
_______________________________________
Brian O'Meara
Assistant Professor
Dept. of Ecology & Evolutionary Biology
U. of Tennessee, Knoxville
http://www.brianomeara.info
Students wanted: Applications due Dec. 15, annually
Postdoc collaborators wanted: Check NIMBioS' website
Calendar: http://www.brianomeara.info/calendars/omeara
On Mon, Sep 17, 2012 at 9:11 PM, Jason S <[email protected]> wrote:
Thanks, guys. That's exactly what I needed.
________________________________
From: Liam J. Revell <[email protected]>
To: Matt Pennell <[email protected]>
-project.org>
Sent: Monday, September 17, 2012 9:22 PM
Subject: Re: [R-sig-phylo] variation in rates over time
Hi Jason. Matt is absolutely correct. You can do this with phytools.
Say, for instance, you have an ultrametric phylogeny with branches in
millions of years (tree) and data vector containing the trait values for
species (x) and you want to test the hypothesis that the last 3.4 my has
a different rate of evolution than the rest of the tree, you could do
this as follows:
library(phytools) # load phytools
tree<-make.era.map(tree,c(0,max(nodeHeights(tree))-3.4))
plotSimmap(tree,pts=F,lwd=3) # visualize
fit<-brownie.lite(tree,x) # fit model
That's it. Good luck. Liam
Liam J. Revell, Assistant Professor of Biology
University of Massachusetts Boston
web: http://faculty.umb.edu/liam.revell/
email: [email protected]
blog: http://phytools.blogspot.com
On 9/17/2012 7:46 PM, Matt Pennell wrote:
Jason,
I think the best way to do this is with the approach of O'Meara et al.
2006
Evolution "Brownie".
Liam Revell has implemented this in R in his package phytools. You can
modify the steps taken in this tutorial here
http://phytools.blogspot.com/2011/07/running-brownielite-for-arbitrarily.html
perhaps
in conjunction with the function make.era.map()
http://phytools.blogspot.com/2011/10/new-function-to-plot-eras-on-tree.html
though
I admittedly have not tried this myself.
Perhaps Liam or someone else has a better explanation but I hope this is
at
least somewhat helpful.
cheers,
matt
Hello,
I see that there are several interesting alternatives to test for
different rates among clades. However, I was wondering if there is a
method
to test for varying rates over time. I'm aware of Pagel's delta and the
EB
model, but I was thinking more in terms of testing if there is a
different
rate for the entire tree after a specified point in time. For instance,
if
a snail predator colonizes an island 3.4 Mya, is there evidence for an
increased rate of evolution in the prey after that point in time?
Something
like two lambdas, one for before and one for after that point in time.
Thanks!
Jason
[[alternative HTML version deleted]]
_______________________________________________
R-sig-phylo mailing list
[email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
[[alternative HTML version deleted]]
_______________________________________________
R-sig-phylo mailing list
[email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
[[alternative HTML version deleted]]
_______________________________________________
R-sig-phylo mailing list
[email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
[[alternative HTML version deleted]]
_______________________________________________
R-sig-phylo mailing list
[email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
_______________________________________________
R-sig-phylo mailing list
[email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
