Hi everyone,
I want to compare Pagel's lambda for phylogenetic signal in a given trait
against randomized lambdas by doing a tip shuffle on my tree of interest and
recalculating lambda each time. I know people often compare an observed lambda
value with a lambda of zero or 1, but I figured this might be another (perhaps
more straightforward?) way to do it. I tried to write a wrapper for the
fitContinuous function in geiger. I would appreciate it if someone could take a
quick look at it and see if it's valid or provide some feedback, as I don't
remember seeing this done before. Thanks so much!
library(geiger)
library(picante)
lambda_shuffle <- function(phylogeny, trait, reps){
observed.lambda <- (fitContinuous(phy=phylogeny, dat=trait,
model="lambda"))$opt$lambda
tipsh <- replicate(reps, (fitContinuous(phy=tipShuffle(phylogeny), dat=trait,
model="lambda"))$opt$lambda)
mean.tipsh <- mean(tipsh)
p.tipsh <- (reps + 2 - rank(c(observed.lambda, tipsh))[1])/(reps+1)
results <- data.frame(observed.lambda, mean.randomized.lambda = mean.tipsh,
P.value = p.tipsh)
return(results)
}
Vincenzo
[[alternative HTML version deleted]]
_______________________________________________
R-sig-phylo mailing list - [email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
Searchable archive at http://www.mail-archive.com/[email protected]/
