Hi Nina- Looking at your data, you have an enormous effect of large values. If you log-transform things, this behaves much better. One obvious flag that suggests a problem is that all of your fitted values are on the "bounds" of fitContinuous. Basically, fitContinuous uses a default method of optimizing parameters that constrains the search to a pre-defined min and max. Both your fitted BM (beta = 20) and lambda (1e-07) were coming up on the boundaries, indicating that the likelihood optimization was failing.
Probably Liam Revell has code in his phytools package (http://phytools.blogspot.com/) that provides more robust fitting than fitContinuous (e.g., that does not rely on pre-defined min-max bounds), or I can send you something if not! But I'm pretty sure you will want to log transform. Also, by testing BM versus lambda as below, you are explicitly testing lambda equal to 1 (pure Brownian motion) versus lambda NOT equal to 1. This is not the same as testing lambda = 0: BM is a model with lambda = 1. I recommend also trying phylosignal in the picante package to help get some intuition for what is going on! ~Dan Rabosky On Nov 12, 2011, at 8:28 AM, Nina Hobbhahn wrote: > Hi there, > > I am using fitContinuous to test whether the number of stomata on orchid > petals exhibits significant phylogenetic signal. Please find my code and data > files below. The stoma count represents species averages. I'm somewhat > confused by the estimate of lambda in the lambdaFit model: It is 1e-07, yet > the when I compare the likelihoods of a lambdaFit model and a brownFit model, > the chi-square test suggests that this lambda value is highly significantly > different from zero. Am I doing something wrong? > > Assuming this result is correct, I presume the correct interpretation would > be that there is significant yet very weak phylogenetic signal? > > Any help would be greatly appreciated. > > Sincerely, > > > Nina Hobbhahn, PhD cand. > > 280 Biological Sciences > 2500 University Drive NW > University of Calgary > Calgary, AB, T2N1N4 > Canada > > phy <-read.tree("allDisa_aged") > > Disa.dat<-read.csv("Stomacount where present1.csv") > #check for missing values in tree and dataset > IN1<-phy$tip.label %in% Disa.dat$species > phy$tip.label[!IN1] > IN2<-Disa.dat$species %in% phy$tip.label > Disa.dat$species[!IN2] > notinphy <- as.vector(Disa.dat$species) > # make new tree > phy2 <- drop.tip(phy, which(!phy$tip.label %in% notinphy)) > > allstomata<-Disa.dat$allstomata > names(allstomata)<-Disa.dat$species > DF.Disa<-data.frame(allstomata,species=Disa.dat$species) > DF.Disa<-DF.Disa[phy2$tip.label, ] > DF.Disa > > # run and store analysis: > brownFit<-fitContinuous(phy2, allstomata) > brownFit > # treetransform lambda > lambdaFit<-fitContinuous(phy2, allstomata, model="lambda") > lambdaFit > #compare likelihoods: > G=2*(lambdaFit$Trait1$lnl-brownFit$Trait1$lnl) > G > P.lambda=pchisq(G,1, lower.tail=FALSE) > P.lambda > > <allDisa_aged.txt> > <Stomacount where present1.xls> > > > > > > _______________________________________________ > R-sig-phylo mailing list > R-sig-phylo@r-project.org > https://stat.ethz.ch/mailman/listinfo/r-sig-phylo [[alternative HTML version deleted]] _______________________________________________ R-sig-phylo mailing list R-sig-phylo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
