Hi Karin.
gls(...,correlation=corPagel(...,fixed=F),method="ML") and
pgls(...,lambda="ML") should produce the same result unless one or the
other is not optimizing properly or unless different bounds on the
fitted model parameters (in particular, lambda) are used.
To show this to yourself, you can try the following example which
simulates and then fits the same model using three different functions.
I hope this is helpful. - Liam
set.seed(1)
require(phytools); require(caper); require(nlme); require(geiger)
# simulate tree
tree<-pbtree(n=100,scale=1)
# simulate x under BM
x<-fastBM(tree)
# simulate y & residual error under the lambda model
y<-0.4*x+fastBM(lambdaTree(tree,0.7),sig2=0.5)
data<-as.data.frame(cbind(x,y))
data<-cbind(names(x),data)
colnames(data)<-c("Species","x","y")
# fit using nlme:gls
r1<-gls(y~x,data=data,correlation=corPagel(0.5,tree,fixed=FALSE),method="ML")
# fit using caper:pgls
r2<-pgls(y~x,data=comparative.data(tree,data,"Species"),lambda="ML")
# fit using phytools:phyl.resid
r3<-phyl.resid(tree,x,y,method="lambda")
--
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 1/13/2012 6:04 AM, Karin Schneeberger wrote:
Dear all
I am currently analysing a set of traits from different species - let's say
brain size and body mass - using a tree where branch lengths were set
artificially (after Graphen) as the length was not available.
I calculate phylogenetic independent contrasts (PIC) using the package "caper"
and the crunch-algorithm. Beforehand, I used the function caic.diagnostic to test whether
the model is appropriate (which was the case). In the end, this provided me with a
significant correlation between brain size and body mass.
I then calculate lambda as a measurement for phylogenetic signal, using the package
"ape":
data<- comparative.data(tree, data, Species)
result<- gls(brainsize~bodymass, data = data, correlation = corPagel(value = 1, phy =
tree, fixed = FALSE), method = "ML")
summary(result)
Here, lambda was 1, which - as far as I understood - means evolution under Brownian
motion. I perform a phylogenetic generalised least square (PGLS) regression using
"caper", which includes lambda, as an addition to the PIC analysis:
data<- comparative.data(tree, data, Species, vcv = TRUE)
Mod1<- pgls(brainsize~bodymass, data = data, lambda = "ML")
summary(Mod1)
Surprisingly, lambda here was zero. Furthermore, the correlation between brain
size and body mass became non-significant.
I now have two contradiction values for lambda and results from both PIC and
PGLS that differ dramatically in the way they can be interpreted. How does it
come that with one method, lambda is zero and with the other lambda is one? And
in the end, which method should be used, PIC or PGLS, or should both be
discussed?
I would be grateful for advice.Thank you very much in advance for your help.
Best regards
Karin Schneeberger
PhD candidate in Ecology and Evolution
[[alternative HTML version deleted]]
_______________________________________________
R-sig-phylo mailing list
R-sig-phylo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
_______________________________________________
R-sig-phylo mailing list
R-sig-phylo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
