Warning using phyloglm

[paternogbc]

Hi,

I am trying to run a logistic phylgenetic regression with 'phyloglm' and would like to understand this warning:

Warning message: In phyloglm(y ~ x, dados, phy = phy) : the boundary of the linear predictor has been reached during the optimization procedure. You can increase this bound by increasing 'btol'.

Even if I try to increase 'btol', the warning persists. Could you please provide some more details about this warning. Should I be concerned about the robustness of the results?

Thanks in advance, Gustavo Paterno

[lamho86]

Dear Gustavo Paterno,

The parameter 'btol' controls the constraint of fitted values between (a,b) where a is close to 0 and b is close to 1. For example, the default btol=10 makes a = 1/(1+exp(10))=0.000045, and b = 1/(1+exp(-10)) = 0.999955.

This warning is to inform that during the optimization search, the boundary has been reached at some points. So, there is a possibility that be optimal solution lies outside of the constraint and the current solution is a point at the boundary which is closest to the optimal solution. That's why phyloglm suggests to increase 'btol'. You could check to see if there is any problem by looking at the fitted values from the output.

Hope it helps, Lam

[paternogbc]

Dear @lamho86,

thanks for your fast reply. I was understanding the warning in that way, but was not so sure. I have increased the Btol but it continues to print this warning, so maybe the optimal solution is realy very close to zero. Actually, what I did was a z-transformation of my X variable and now there is no more warnings. Do you see any problem with this approach?

I will give a close look in the predicted values to check if there is any problem.

Thanks once again!

All the best, Gustavo

[lamho86]

Hi Gustavo,

The warning does not mean that the optimal solution is really close to either 0 or 1. It just means that there is that possibility. So, you have to check.

Regarding the transformation, I don't see any problem as long as you can still interpret your results.

Best, Lam

[paternogbc]

Hi Lam,

thanks for your explanation.

We have developed an R package for sensitivity analysis in comparative methods sensiPhy which is mainly based on phylolm functions. Specially because you code is much faster then other packages. But also because phylolm functions allow for other types of evolutionary models when fitting PGLS. We are going to release a stable version to CRAN on dezember this year.

All the best, Guga

[paternogbc]

Dear Lam,

If you dont mind, I would like to ask a second question. How can I extract the phylogenetic residuals from models fitted with phylolm. As I understand, this function only provides the RAW residuals.

Do you have any advice on how can I get this values using phylolm package?

All the best, Gustavo

[cecileane]

Congratulation on sensiPhy Gustavo!

The residuals returned by phylolm are raw, i.e. they are the differences between the observed values and the expected values. They are still phylogenetically correlated. Do you mean that you want phylogenetically uncorrelated residuals? For this you need to multiply the raw residuals by the inverse square root of the covariance matrix. There is a fast way to get this matrix, using just one tree traversal, available in this package: https://github.com/khabbazian/l1ou/ see function sqrt_OU_covariance(). So far, this function is restricted to a BM model, or to an OU model on an ultrametric tree. Cecile.

[paternogbc]

Dear @cecileane, thanks for that. sensiPhy is a work in progress.. but we should finish very soon and submit the first version to CRAN in dez. But thank you for this great package (phylolm) it is really very fast and more flexible when compared to any other package available!

Regarding the residuals, thanks for clarifying. Yep, what I want is the uncorrelated phylogenetic residuals.

Actually, I want to test if my model residuals still has phylogenetic signal after fitting the PGLS. I am not sure what residuals I should use to test that. Thanks for suggesting the covariance matrix approach, I will try that. Could you confirm if this approach would give the same output as the mod$phyres from the pgls function in caper (like the code below):

mod <- pgls(formula, data) mod$phyres

Thanks again,

all the best, Guga

[cecileane]

There are many ways to uncorrelate the residuals: all of them use some inverse square root C=V^{-1/2} of the covariance V. There are many square-roots of V, actually: there are many matrices C such that CVC' = I, so there are many ways to uncorrelate the residuals. There is no unique way.

One standard choice is to constrain C to be symmetric, and with positive eigenvalues. Perhaps that's what coded in the pgls function in caper, I have not looked at the details.

Another choice here, and that's the choice we made in the l1ou package, is a matrix V^{-1/2} (also a matrix V^{+1/2}) such that each row in these matrices corresponds to an internal node of the tree. That gives more ways to test for model adequacy: to check that the phylogenetically corrected residuals do not actually correlate with node height, or other node attribute.

Anyway, the short version is that two sets of phylogenetically-corrected residuals need not be the same.

[paternogbc]

Dear @cecileane, thank you so much for this explanation!

All the best, Guga

[lamho86]

Hi Gustavo,

If you do a small number of tests, you can compute this C=V^{-1/2} matrix using any standard matrix decomposition function in R (for example svd()) because you only compute it once for each test.

Lam

On Tue, Oct 20, 2015 at 11:16 AM, Gustavo Paterno notifications@github.com wrote:

Dear @cecileane https://github.com/cecileane, thank you so much for this explanation!

All the best, Guga

— Reply to this email directly or view it on GitHub https://github.com/lamho86/phylolm/issues/3#issuecomment-149652604.

[paternogbc]

Hi @lamho86 ,

thanks for your suggestion. I will have a look.

all the best, Guga