eliotmiller / metricTester

8 stars 3 forks source link

expectations() fails for metric "PD_Cadotte" #5

Closed sagitaninta closed 7 years ago

sagitaninta commented 7 years ago

Dear Dr. Miller,

I am using metricTester to explore null model behavior of my data set for my research and encountering errors from using expectations() function. I am interested in the regional null models proposed in the Ecography paper and wanted to know how my data set behave in this null model using MPD and PD metrics but it gives error every time I choose "PD_Cadotte" for metrics argument. Here is my minimum reproducible example using the code in the example.

tree <- geiger::sim.bdtree(b=0.1, d=0, stop="taxa", n=50)
sim.abundances <- round(rlnorm(5000, meanlog=2, sdlog=1)) + 1
cdm <- simulateComm(tree, richness.vector=10:13, abundances=sim.abundances)

test <- expectations(picante.cdm=cdm, tree=tree, optional.dists=NULL,
regional.abundance=NULL, distances.among=NULL, randomizations=3, cores="seq",
nulls="old_regional", metrics=c("PD_Cadotte", "NAW_MPD"),
concat.by="both", output.raw=FALSE)
Error in data.frame(richness = null.output$richness, metric = null.output[,  : 
  arguments imply differing number of rows: 0, 12
In addition: Warning messages:
1: Not running analysis in parallel. See 'cores' argument. 
2: In prepNulls(tree, picante.cdm, regional.abundance, distances.among) :
  Regional abundance not provided. Assumed to be equivalent to CDM
3: In prepNulls(tree, picante.cdm, regional.abundance, distances.among) :
  Distances among plots not provided. Null models that require this input will not be run

The example code run as expected (nulls="richness", metrics=c("richness", "NAW_MPD")), but every time I use "PD_Cadotte" for any null models, the same error message occurs. Otherwise, it works for metrics and null interests of my choice... so far. I cannot use "PD" in the metrics argument as my tree is not rooted.

Does this mean that the choices in argument nulls and metrics cannot be deliberately used for any list of metrics and null models in the package?

Thank you in advance for your response.

eliotmiller commented 7 years ago

Hi @sagitaninta. PD_Cadotte is just a slightly modified form of phylogenetic diversity that includes the root of the tree. metricTester accesses it via the picante function pd(), with the include.root argument set to TRUE. If your tree is unrooted, you can't calculate that form of PD because it doesn't "make sense". That said, I think you could still calculate the original form of PD, if you think it's biologically meaningful for your case. In practice, these two forms of PD are usually approximately the same thing. Here's an example of what I'm talking about

library(picante)
data(phylocom)
# This works whether you include the root or not, because the tree is rooted
pd(phylocom$sample, phylocom$phylo, include.root=TRUE)
pd(phylocom$sample, phylocom$phylo, include.root=FALSE)
# This only works when the tree is rooted, because you obviously can't
# include the root in the calculation if there isn't one
pd(phylocom$sample, unroot(phylocom$phylo), include.root=FALSE)
pd(phylocom$sample, unroot(phylocom$phylo), include.root=TRUE)

Does that help?

eliotmiller commented 7 years ago

@sagitaninta An additional thing I just noticed. I say it in the documentation somewhere, but I need to improve the way metricTester is coded so that you can't make this mistake--it's easy to miss. Any time you want to concatenate your randomized results by species richness, you need to have included that as a metric. In your example above, you ask to summarize the randomizations both by richness and by plot (concat.by = "both"). But, then you don't tell it to calculate species richness as one of the metrics. Again, I totally admit this isn't intuitive and I'll try to find a minute to make it easier. In the meantime, if you change your example above to include "richness" as a metric, it works just fine:

tree <- geiger::sim.bdtree(b=0.1, d=0, stop="taxa", n=50)
sim.abundances <- round(rlnorm(5000, meanlog=2, sdlog=1)) + 1
cdm <- simulateComm(tree, richness.vector=10:13, abundances=sim.abundances)

test <- expectations(picante.cdm=cdm, tree=tree, optional.dists=NULL,
regional.abundance=NULL, distances.among=NULL, randomizations=3, cores="seq",
nulls="old_regional", metrics=c("richness", "PD_Cadotte", "NAW_MPD"),
concat.by="both", output.raw=FALSE)