Hi Megan, This looks to me like the problems are arising because the optimal value of lambda is negative (and you are running into optimization problems by bounding it at zero). This is a bit strange -- but not unheard of, I have seen similar results before. What it implies that the trait values of closely related species are more different than one would suspect from random. One might imagine that a process like character displacement could generate such a pattern (though this is not evidence that character displacement is necessarily important).
It is difficult to know what to do here as comparative methods haven't really been developed to deal with this case. I don't think it makes sense to not set bounds on lambda as I am not sure how to interpret a negative branch length transformation (which is what a negative lambda means). So as gls will equal pgls if lambda=0, I would just run regular regression without correcting for phylogeny (heresy, I know!). I would also look at the residuals of individual data points to look for outliers which may be driving this weird result (i.e. if two very closely related species have very different trait values). A simple way to test this would be to compute the independent contrasts (using pic()) and then see if there are any values that are way bigger than the rest. Not sure what else to tell you here. Best of luck with your analysis. Hope this was at least somewhat helpful. matt On Tue, Jun 11, 2013 at 12:25 PM, Megan Bartlett <mkba...@gmail.com> wrote: > My apologies, I left off the profile graph. Attached here! > > Thanks very much. > > > On Tue, Jun 11, 2013 at 12:19 PM, Megan Bartlett <mkba...@gmail.com>wrote: > >> Hi all, >> >> I'm replying to an old thread ( >> https://stat.ethz.ch/pipermail/r-sig-phylo/2012-January/001813.html), >> because I'm having a very similar problem but I wasn't quite sure from the >> answers how to resolve it. I'm trying to fit a PGLS model between trait >> variables using both Lambda and OU models for phylogeny, and to determine >> whether those models are better fits than a phylogenetically-uncorrected >> linear model using AICc values. I was originally using gls() and REML >> optimization, which returned seemingly reasonable estimates for lambda. I >> then realized that REML isn't suitable for AICc comparison and specified >> method = "ML", and I can't calculate results for lambda for some (but not >> all) traits- I get the following error message no matter what I specify for >> initial lambda values in corPagel(): >> >> Error in gls(TLP_DS ~ AverageOP_DS, PhyloDataUN, pagelcor, method = "ML") >> : >> false convergence (8) >> >> I then tried pgls() in caper, and I get a similar shape for >> pgls.profile() for every one of my pairwise trait correlations (profile >> attached). Almost every best-fit lambda is the lowest boundary for lambda >> estimation, with a few exceptions. Every PGLS model produces a lambda value >> where the lower 95% CI can't be calculated (returns NA), probably because >> log likelihood either increases or remains constant as the lambda values >> used approach 0. Does that mean that >> >> 1) The pgls() ML search is failing like the gls() optimization, and >> that's why it's returning lower 95%CI = NA and lambda = 1E-6? >> >> or 2) The pgls estimates are reasonable and I can interpret the lower >> 95%Ci = NA result as none of the lambda values being significantly greater >> than 0, thus phylogenetic history does not significantly affect my >> univariate trait correlations. >> >> and 3) If pgls() is failing despite what I specify for initial >> parameters, is there anything else I can do to correct for phylogeny? >> >> Thanks so much for your help! I'll post example code below (please let me >> know if I can supply sample data, I didn't want to lengthen an already long >> post). >> >> Best, >> >> Megan Bartlett >> >> ## PGLS model that returned a result >> >> comparative.data(XTBGTree, PhyloDataUN, names.col = "Species", vcv =TRUE, >> vcv.dim=3, warn.dropped=TRUE) -> ComparedData >> >> pgls(TLP_DS ~ AverageOP_DS, ComparedData, lambda = "ML") -> model1 >> >> pgls.profile(model1, "lambda", N =100, param.CI = 0.95) -> modelprofile >> >> >> #Call: >> >> #pgls(formula = TLP_DS ~ AverageOP_DS, data = ComparedData, lambda = "ML") >> >> >> #Residuals: >> >> # Min 1Q Median 3Q Max >> >> #-0.0021551 -0.0004180 -0.0000137 0.0002800 0.0033862 >> >> >> #Branch length transformations: >> >> >> #kappa [Fix] : 1.000 >> >> #lambda [ ML] : 0.000 >> >> #lower bound : 0.000, p = 1 >> >> #upper bound : 1.000, p = 9.0506e-05 >> >> #95.0% CI : (NA, 0.418) >> >> #delta [Fix] : 1.000 >> >> >> #Coefficients: >> >> #Estimate Std. Error t value Pr(>|t|) >> >> #(Intercept) 0.151049 0.019080 7.9168 8.64e-10 *** >> >> 3AverageOP_DS 1.751522 0.017935 97.6569 < 2.2e-16 *** >> >> #--- >> >> #Signif. codes: 0 *** 0.001 ** 0.01 * 0.05 . 0.1 1 >> >> >> #Residual standard error: 0.0009188 on 41 degrees of freedom >> >> #Multiple R-squared: 0.9957, Adjusted R-squared: 0.9956 >> >> #F-statistic: 9537 on 2 and 41 DF, p-value: < 2.2e-16 >> >> >> ###GLS version that did not: >> >> corPagel(0.1, XTBGTree, fixed =FALSE)-> pagelcor >> >> gls(TLP_DS ~ AverageOP_DS, PhyloDataUN, pagelcor, method = "ML") -> >> pagel.model >> Error in gls(TLP_DS ~ AverageOP_DS, PhyloDataUN, pagelcor, method = "ML") >> : >> false convergence (8) >> >> > > _______________________________________________ > 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/ > [[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/
