Apparently the plot did not come through correctly. I'm attaching a PDF of it.

Best, Brian _______________________________________________________________________ Brian O'Meara, http://www.brianomeara.info, especially Calendar <http://brianomeara.info/calendars/omeara/>, CV <http://brianomeara.info/cv/>, and Feedback <http://brianomeara.info/teaching/feedback/> Associate Professor, Dept. of Ecology & Evolutionary Biology, UT Knoxville Associate Head, Dept. of Ecology & Evolutionary Biology, UT Knoxville Associate Director for Postdoctoral Activities, National Institute for Mathematical & Biological Synthesis <http://www.nimbios.org> (NIMBioS) On Wed, Apr 4, 2018 at 4:24 PM, Brian O'Meara <omeara.br...@gmail.com> wrote: > Hi, Rafael et al. > > ?OUwie has some information on the standard errors: > > The Hessian matrix is used as a means to estimate the approximate standard >> errors of the model parameters and to assess whether they are the maximum >> likelihood estimates. The variance-covariance matrix of the estimated >> values of alpha and sigma^2 are computed as the inverse of the Hessian >> matrix and the standard errors are the square roots of the diagonals of >> this matrix. The Hessian is a matrix of second-order derivatives and is >> approximated in the R package numDeriv. So, if changes in the value of a >> parameter results in sharp changes in the slope around the maximum of the >> log-likelihood function, the second-order derivative will be large, the >> standard error will be small, and the parameter estimate is considered >> stable. On the other hand, if the second-order derivative is nearly zero, >> then the change in the slope around the maximum is also nearly zero, >> indicating that the parameter value can be moved in any direction without >> greatly affecting the log-likelihood. In such situations, the standard >> error of the parameter will be large. >> >> For models that allow alpha and sigma^2 to vary (i.e., OUMV, OUMA, and >> OUMVA), the complexity of the model can often times be greater than the >> information that is contained within the data. As a result one or many >> parameters are poorly estimated, which can cause the function to return a >> log-likelihood that is suboptimal. This has great potential for poor model >> choice and incorrect biological interpretations. An eigendecomposition of >> the Hessian can provide an indication of whether the search returned the >> maximum likelihood estimates. If all the eigenvalues of the Hessian are >> positive, then the Hessian is positive definite, and all parameter >> estimates are considered reliable. However, if there are both positive and >> negative eigenvalues, then the objective function is at a saddlepoint and >> one or several parameters cannot be estimated adequately. One solution is >> to just fit a simpler model. Another is to actually identify the offending >> parameters. This can be done through the examination of the eigenvectors. >> The row order corresponds to the entries in index.matrix, the columns >> correspond to the order of values in eigval, and the larger the value of >> the row entry the greater the association between the corresponding >> parameter and the eigenvalue. Thus, the largest values in the columns >> associated with negative eigenvalues are the parameters that are causing >> the objective function to be at a saddlepoint. >> > > You can get better (by which I mean more accurate) estimates by > parametric bootstrapping (simulation). > > Doing a quick summary of your OUM results (code at end): > > Trait: contr > nonfor: 0.12 (-0.12, 0.36) > inter: 0.12 (-0.12, 0.36) > fores: 0.13 (-0.13, 0.39) > Phylogenetic half life: 0.06 > Tree height: 28.34 > > Trait: cshd > nonfor: 0.24 (-0.23, 0.72) > inter: 0.33 (-0.32, 0.98) > fores: 0.2 (-0.19, 0.6) > Phylogenetic half life: 6.35 > Tree height: 28.34 > > Trait: dors > nonfor: 0.04 (-0.04, 0.12) > inter: 0.15 (-0.14, 0.44) > fores: 0 (0, 0.01) > Phylogenetic half life: 3.68 > Tree height: 28.34 > > Trait: vent > nonfor: 0.05 (-0.05, 0.14) > inter: 0.32 (-0.3, 0.93) > fores: -0.03 (0.03, -0.08) > Phylogenetic half life: 5.43 > Tree height: 28.34 > > It does seem that you don't have a lot of ability to tell optima apart. > The alpha value is pretty low, though not tremendously (the half life is a > good way to look at that). But the raw data also suggest it's not three > very clear clumps (beeswarm plot, regimes on the x axis, trait 1 value on > the y). It's not phylogenetically corrected, but still a good gut check: > > > > I agree that for your data, it may be that simpler models work better. It > could also explain why you're getting crazy theta values for models with V > or A allowed to vary: without strong attraction, it's hard to figure out > what the attraction is towards, and especially for the third regime, with > only 11 data points, there probably isn't a lot of info to estimate many > parameters for just that regime. > > Best, > Brian > > library(OUwie) > library(ape) > library(beeswarm) > > data <- read.csv("~/Downloads/OUM_PARAMS_ant_f.csv", row.names=1) > load("~/Downloads/examplemodel.Rdata") > > traits <- colnames(data) > > regimes <- c("nonfor", "inter", "fores") > > observations <- list() > > for (trait.index in sequence(length(traits))) { > cat(paste("\n\nTrait:", traits[trait.index])) > for (regime.index in sequence(length(regimes))) { > point.estimate <- data[paste0("theta.",regimes[regime.index]), > traits[trait.index]] > se <- data[paste0("theta.",regimes[regime.index]), > traits[trait.index]] > result.vector <- round(c(point.estimate,point.estimate-1.96*se, > point.estimate+1.96*se),2) > cat(paste0("\n",regimes[regime.index], ": ", result.vector[1], " > (", result.vector[2], ", ", result.vector[3], ")")) > if (trait.index==1) { > observations[[regime.index]] <- > bla$data[which(bla$data[,1]==regime.index),2] > > } > } > cat(paste0("\nPhylogenetic half life: ", round(log(2)/data["alpha", > traits[trait.index]], 2))) > cat(paste0("\nTree height: ", round(max(ape::branching. > times(bla$phy)),2))) > } > > points.to.plot <- bla$data[,1:2] > colnames(points.to.plot) <- c("regime", "trait1") > beeswarm(trait1~regime, data=points.to.plot, col="black", pch=20) > > > > > _______________________________________________________________________ > Brian O'Meara, http://www.brianomeara.info, especially Calendar > <http://brianomeara.info/calendars/omeara/>, CV > <http://brianomeara.info/cv/>, and Feedback > <http://brianomeara.info/teaching/feedback/> > > Associate Professor, Dept. of Ecology & Evolutionary Biology, UT Knoxville > Associate Head, Dept. of Ecology & Evolutionary Biology, UT Knoxville > Associate Director for Postdoctoral Activities, National Institute for > Mathematical & Biological Synthesis <http://www.nimbios.org> (NIMBioS) > > > On Wed, Apr 4, 2018 at 4:07 PM, Jacob Berv <jakeberv.r.sig.ph...@gmail.com > > wrote: > >> Even if those models may fit the data better, they may not necessarily >> help Rafael determine whether or not the parameter estimates of interest >> are different across regimes (though perhaps BMS might be informative). >> >> Rafael, maybe you could try fixing the ancestral regimes to match most >> likely states for each node from your SIMMAP summary? I wonder if that >> might decrease your ‘uncertainty’ in parameter estimates. >> >> I don’t have a great answer for your main question though, which is a >> good one. >> >> Jake >> >> > On Apr 4, 2018, at 8:59 PM, William Gearty <wgea...@stanford.edu> >> wrote: >> > >> > Hi Rafael, >> > >> > Have you tried running the BM1, BMS, or OU1 models? If so, how do the >> AICc >> > values for those compare to those of the OUM model? If not, you should >> make >> > sure you run those. >> > If you find that the these models fit your data better, this could >> indicate >> > that there isn't a large different across regimes and an OUM model isn't >> > necessarily suitable. >> > >> > Best, >> > Will >> > >> > On Wed, Apr 4, 2018 at 12:30 PM, Rafael S Marcondes < >> raf.marcon...@gmail.com >> >> wrote: >> > >> >> Dear all, >> >> >> >> I'm writing (again!) to ask for help interpreting standard errors of >> >> parameter estimates in OUwie models. >> >> >> >> I'm using OUwie to examine how the evolution of bird plumage color >> varies >> >> across habitat types (my selective regimes) in a tree of 229 tips. I >> was >> >> hoping to be able to make inferences based on OUMV and OUMVA models, >> but I >> >> was getting nonsensical theta estimates from those. So I've basically >> given >> >> up on them for now. >> >> >> >> But even looking at theta estimates from OUM models, I'm getting really >> >> large standard errors, often overlapping the estimates from other >> selective >> >> regimes. So I was wondering what that means exactly. How are these >> erros >> >> calculated? How much do high errors it limit the biological inferences >> I >> >> can make? I'm more interested in the relative thetas across regimes >> than on >> >> the exact values (testing the prediction that birds in darker habitats >> tend >> >> to adapt to darker plumages than birds in more illuminated habitats). >> >> >> >> I have attached a table averaging parameter estimates and errors from >> >> models fitted across a posterior distribution of 100 simmaps for four >> >> traits; and one exemplar fitted model from one trait in one of those >> >> simmaps. >> >> >> >> Thanks a lot for any help, >> >> >> >> >> >> *--* >> >> *Rafael Sobral Marcondes* >> >> PhD Candidate (Systematics, Ecology and Evolution/Ornithology) >> >> >> >> Museum of Natural Science <http://sites01.lsu.edu/wp/mns/> >> >> Louisiana State University >> >> 119 Foster Hall >> >> Baton Rouge, LA 70803, USA >> >> >> >> Twitter: @brown_birds <https://twitter.com/brown_birds> >> >> >> >> >> >> _______________________________________________ >> >> 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-ph...@r-project.org/ >> >> >> >> >> > >> > >> > -- >> > William Gearty >> > PhD Candidate, Paleobiology >> > Department of Geological Sciences >> > Stanford School of Earth, Energy & Environmental Sciences >> > williamgearty.com >> > >> > [[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-ph...@r-project.org/ >> >> _______________________________________________ >> 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-ph...@r-project.org/ >> > >

**
swarm.pdf**

*Description:* Adobe PDF document

_______________________________________________ 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/