Hi All, I am running a pgls regression analysis in the caper package in R. My dataset has a continuous response variable and five explanatory variables (three continuous, two discrete and with multiple categories). Each datapoint in the response variables corresponds to a single species of which there are 67 in total. I plan to run multiple regression models, each model containing a different combination of the explanatory variables (fixed effects only, no interactions; total number of models in full set = 32), and then use AICc values to calculate AICc model weights as per Burnham & Anderson (2002). These weights will then be used to calculate weighted-average regression coefficients for each explanatory variable across all models (also as per Burnham & Anderson 2002). The weighted-average regression coefficients will then be used for prediction. I am not necessarily interested in the extent of the phylogenetic effect - the aim is to compare results from the PGLS models to the same! set of models run using simple OLS regression and (hopefully) show that the results are broadly similar. I have a phylogeny (phyl) of the 67 species taken from the literature. This phylogeny has no branch lengths, so I have set all branch lengths as equal to 1 using: phylC <- compute.brlen(phyl, 1) My questions are: Firstly, should I be optimising the value of the lambda transformation in each of the 32 models to control for the different levels of phylogenetic effect when including different combinations of explanatory variables? e.g. model <- pgls(y ~ x1 + x2 + x3, data=phylC, lambda = "ML") Or would it be sufficient/robust to set lambda = 1 (the default) for all models, thereby simply assuming that there is a phylogenetic effect without necessarily being interested in it's extent. Optimising lambda for all the 32 models returns maximum likelihood values of lambda that range from 0.478 (significantly > 0 and < 1 using likelihood ratio tests) to 0. Secondly, if I should be optimising lambda in each model, is it meaningful to then compare and average across models using AICc values? Surely I could not be sure whether the difference in AICc (i.e. model fit) between models was due to the inclusion of different explanatory variables (which I am really interested in) or the different value of the lambda transformation? In a previous post, David Orme wrote: "...it isn't necessarily safe to use AIC or deletion tests to compare models with different covariance structures _and_ different fixed effects, since changing the branch lengths affects the AIC of the model...Basically, it isn't easy to do model simplification of your fixed effects whilst also finding a ML evolutionary model for that changing set of fixed effects. Another argument for keeping the branch length transforms simple." http://www.mail-archive.com/r-sig-phylo@r-project.org/msg02060.html Though this implies that optimising lambda and then comparing models with AIC is not a good idea, I have found various references where this appears to have been done, so am confused as to what to do. Many thanks in advance for your help, Lucy Malpas Conservation Scientist Species Monitoring and Research, Conservation Science The Royal Society for the Protection of Birds [[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/
