Hi,
I have some questions about the method from this paper:
Paradis et al. (2013) Quantifying variation in speciation and extinction
rates with clade data. DOI: 10.1111/evo.12256
#Some hypothetical data:
SpeciesCount <- c(2,13,2,2,8,1,16)
Ages <- c(1,2,0.5,0.6,2,0.2,3)
#And using the simplest model, assuming extinction rate = 0:
halfdev <- function(p) {
if (p <= 0 || p >= 2) return(1e100)
-sum(dyule(SpeciesCount, p, Ages, log = TRUE))
}
#optimizing with nlm:
res <- nlm(halfdev, .1, hessian = TRUE)
#Questions:
#The halfdev function is a -lnL function, so using the nlm function will
maximize
#the likelihood, correct?
#So, res$estimate is the estimate of lambda (for all clades) with the
highest likelihood?
#This means that I can use the following approach to get clade specific
estimates of speciation rate:
lambdas <- -dyule(SpeciesCount, res$estimate, Ages, log = TRUE)
# lambdas summarized with data
cbind(lambdas,SpeciesCount,Ages)
# lambdas SpeciesCount Ages
# 1.4798877 2 1.0
# 3.6698954 13 2.0
# 1.4207890 2 0.5
# 1.3906988 2 0.6
# 3.0143270 8 2.0
# 0.2096531 1 0.2
# 3.8052613 16 3.0
#And I can use these lambda values as a response variable in a pgls (or
#equivalent) for downstream analysis of speciation rates?
Thanks in advance for all/any help.
Jostein
[[alternative HTML version deleted]]
_______________________________________________
R-sig-phylo mailing list - R-sig-phylo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
Searchable archive at http://www.mail-archive.com/r-sig-phylo@r-project.org/
