Invalid extinction probability, set to 1

[cintiaiha]

Dear Dr Zwaenepoel

I am trying to perform WGD inference using reversible-jump MCMC using the gene count file from OrthoFinder. I am using julia-1.5.2 in a Ubuntu (18.04) remote server.
I tested with the example files 9dicots and the script worked perfectly. However, when I try to run my own data, the logpdf! output is -Inf and I end up with these warnings:

┌ Warning: Invalid extinction probability, set to 1
│ p = NaN
└ @ Beluga ~/.julia/packages/Beluga/IeAJt/src/model.jl:445
┌ Warning: Invalid extinction probability, set to 1
│ p = NaN
└ @ Beluga ~/.julia/packages/Beluga/IeAJt/src/model.jl:445
┌ Warning: Invalid extinction probability, set to 1
│ p = NaN
└ @ Beluga ~/.julia/packages/Beluga/IeAJt/src/model.jl:445
┌ Warning: Invalid extinction probability, set to 1
│ p = NaN
└ @ Beluga ~/.julia/packages/Beluga/IeAJt/src/model.jl:445
┌ Warning: Invalid extinction probability, set to 1
│ p = NaN
└ @ Beluga ~/.julia/packages/Beluga/IeAJt/src/model.jl:445
┌ Warning: Invalid extinction probability, set to 1
│ p = NaN
└ @ Beluga ~/.julia/packages/Beluga/IeAJt/src/model.jl:445
┌ Warning: Invalid extinction probability, set to 1
│ p = NaN
└ @ Beluga ~/.julia/packages/Beluga/IeAJt/src/model.jl:445

I also did the trick to enable NaN safe but this problem persists. I tried the same script a couple of months ago and I did not have this problem.

I was wondering if you know if this is a problem with my data or not.

Thank you
Cintia

[arzwa]

Hi, this might have something to do with priors, initial parameters and the time scale of your phylogeny. In the nine dicots example the tree has branch lengths on a scale of 100 My, and the priors etc. are defined with respect to that, so make sure your prior settings match reasonable values for the time scale of your phylogeny.

If that would not be the issue, I'd suggest to have a try with a subset of the data and see if the problem persists. Also did you filter the gene families so that there is at least one gene in each clade stemming from the root of the species tree?

[fayweili]

Hi! I've encountered a similar problem. My root is at ~280MYA. Any tip on adjusting the priors? Thank you!!

[arzwa]

I realize I was probably to cryptic (and wrong it seems) in my above comment. What I meant is that in the 9dicots example in the repository, the unit of the branch lengths is 1 Gy, i.e.

((vvi:0.11,((ath:0.0682,cpa:0.0682):0.0378,(mtr:0.1,ptr:0.1):0.006):0.004):0.00703,((bvu:0.0433,cqu:0.0433):0.0682,(ugi:0.0835,sly:0.0835):0.028):0.0055);

so the Vitis (vvi) branch length is 110 My here in this case, represented as 0.11. The reason is that for branch lengths on this scale, the rate parameters tend to be in the range [0.01 - 10] or so.