Hi Arne,
Just calculating the slope is straightforward. For tree and column
vectors x & y (in order tree$tip.label):
> C<-vcv.phylo(tree)
> ax<-sum(solve(C,x))/sum(solve(C))
> ay<-sum(solve(C,y))/sum(solve(C))
> beta1<-sqrt(t(y-ay)%*%solve(C,y-ay)/(t(x-ax)%*%solve(C,x-ax)))
The model intercept can be obtained by plotting the slope through the
phylogenetic mean of x & y.
> beta0<-ay-beta1*ax
I wrote a function to do this; but also to (optionally) fit Pagel's
lambda jointly for x & y and to return the residuals in y. I did this
last year, and then fixed a few bugs/errors when I saw your query this
morning. I have put it online here:
http://anolis.oeb.harvard.edu/~liam/R-phylogenetics/phy.RMA/v0.1/phyl.RMA.R.
Note, that this *does not* test hypotheses about the RMA regression nor
does it incorporate individual variation (as in the Ives et al. 2007
article). Adding the former would be fairly simple using equations 15
of Ives et al., I think (I will do this when I get a chance); adding the
latter not so much. I believe Tony Ives may have MATLAB code that does
this already. Maybe he can comment.
I hope this is somewhat helpful.
Sincerely, 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 4/20/2011 2:27 AM, Arne Mooers wrote:
Dear members:
Does anyone know of scripts to both estimate and test
phylogenetically-corrected RMA regression slopes (perhaps using the relevant
equations from Ives et al. (Syst. Biol. 2007))?
Cheers,
Arne Mooers
______________________________
Dr. Arne Mooers
Biological Sciences,
Simon Fraser University
Burnaby, BC Canada V5A 1S6
tel. +1 778 782 3979
skype: arnemooers
http://www.sfu.ca/~amooers
http://www.sfu.ca/fabstar
http://www.scientists-4-species.org
_______________________________________________
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
