Huge differences in log marginal likelihoods

[weedy23]

Hello,

I am doing model selection using the generalized stepping stone method to compare the use of a strict vs random local clock models. I have struggled to get ESS > 200 for some posterior estimates of some parameters (some divergence time estimates and the rateChangeCount for the random local clock) with the random local clock model so I predict that the strict clock will be preferred.

Having conducted the analysis, indeed the strict clock is (apparently overwhelmingly) selected as the best model. However the log marginal likelihood value for the random clock model is huge (-3.44 x 10e301) compared to that of the strict clock (-11373). So I am wondering whether this result can be trusted or whether something is not right. 

Can anyone suggest what might be going on here and whether this result can be trusted? Below is the output from the random local clock MLE test.  

Many thanks in advance

Sara

PathParameter        MaxPathLikelihood        MLContribution
0.0000        -1.5995e+308        -3.4460e+301
2.1544e-07        -2.4257e+08        -481.39
2.1715e-06        -2.4735e+07        -160.71
8.3895e-06        -6.3594e+06        -92.619
2.1888e-05        -2.4002e+06        -64.809
4.6050e-05        -1.1813e+06        -52.024
8.4561e-05        -6.6031e+05        -43.595
0.00014136        -3.7861e+05        -35.583
0.00022061        -2.6370e+05        -33.409
0.00032669        -1.7560e+05        -29.958
0.00046416        -1.3458e+05        -28.040
0.00063774        -93118        -26.046
0.00085232        -79727        -25.710
0.0011130        -66423        -24.846
0.0014248        -52181        -24.009
0.0017932        -45511        -23.592
0.0022237        -37968        -23.243
0.0027216        -33842        -23.206
0.0032929        -30592        -23.440
0.0039432        -26502        -23.333
0.0046784        -23895        -23.412
0.0055047        -22732        -24.250
0.0064280        -21575        -24.817
0.0074546        -19859        -25.849
0.0085909        -18546        -26.431
0.0098431        -17648        -27.324
0.011218        -16997        -28.482
0.012722        -16405        -29.826
0.014361        -16064        -30.963
0.016143        -15461        -32.313
0.018075        -14893        -34.034
0.020162        -14583        -35.407
0.022413        -14325        -36.977
0.024834        -14084        -38.916
0.027432        -13748        -40.850
0.030215        -13611        -42.639
0.033190        -13362        -44.846
0.036364        -13146        -46.823
0.039745        -13162        -49.117
0.043339        -12821        -51.270
0.047156        -12801        -53.871
0.051201        -12738        -56.334
0.055484        -12665        -59.084
0.060011        -12504        -61.618
0.064790        -12438        -64.600
0.069830        -12325        -67.303
0.075138        -12281        -70.178
0.080722        -12246        -73.345
0.086591        -12155        -76.515
0.092751        -12105        -79.766
0.099213        -11998        -83.204
0.10598        -11991        -86.592
0.11307        -11962        -90.202
0.12048        -11877        -93.841
0.12823        -11881        -97.643
0.13631        -11881        -101.60
0.14475        -11803        -105.45
0.15355        -11787        -109.48
0.16271        -11773        -113.69
0.17226        -11717        -117.82
0.18218        -11710        -122.34
0.19250        -11686        -126.71
0.20322        -11684        -131.28
0.21436        -11607        -135.83
0.22591        -11639        -140.68
0.23789        -11615        -145.53
0.25031        -11567        -150.46
0.26318        -11559        -155.50
0.27650        -11557        -160.65
0.29029        -11535        -165.96
0.30455        -11530        -171.26
0.31930        -11494        -176.67
0.33454        -11482        -182.15
0.35028        -11478        -187.83
0.36653        -11453        -193.44
0.38330        -11416        -199.04
0.40060        -11407        -204.84
0.41844        -11395        -210.88
0.43683        -11392        -217.38
0.45578        -11408        -223.99
0.47530        -11369        -229.93
0.49539        -11341        -236.06
0.51608        -11344        -242.70
0.53735        -11324        -249.39
0.55924        -11315        -256.10
0.58174        -11304        -262.92
0.60487        -11293        -269.83
0.62863        -11288        -276.94
0.65304        -11286        -284.11
0.67811        -11271        -291.52
0.70384        -11269        -298.90
0.73025        -11264        -306.44
0.75734        -11244        -313.95
0.78513        -11253        -321.85
0.81363        -11234        -329.57
0.84284        -11223        -337.42
0.87278        -11213        -345.40
0.90345        -11212        -353.73
0.93488        -11193        -361.80
0.96705        -11191        -370.31
1.0000        -11194        

log marginal likelihood (using generalized stepping stone sampling) from (pathLikelihood.source - pathLikelihood.destination) = -3.445973394038929E301

[Yasmin97]

Hi Sara, 

That is a pretty huge difference in marginal likelihoods. If I were you, I would double check my data sets or the analysis settings? I have never seen such a massive difference between models! At most, I have seen a difference of perhaps a few hundred between two models (however, I have never compared a strict to a random clock). 

Do you have the output from the generalised-stepping-stone analysis? 

Best,

Yasmin 

[David Swofford]

This result *cannot* be trusted: -10^301 comes very close to the largest (negative) number that can be represented in floating-point format on a 64-bit machine (~ -10^308). A negative log likelihood of that magnitude implies a likelihood of (effectively) zero. A likelihood near zero may indicate a problem with the priors in some part(s) of the parameter space. If the chain wanders into region where the priors have a probability mass near zero, an MCMC sampler will go haywire and will probably fail to converge (ever).

So something has definitely gone wrong, either with your input data/commands, an issue or bug in the MCMC sampling, or an invalid choice of priors. Unfortunately, I'll have to leave it to actual BEAST experts to help you figure out what's going on.

Dave

On Feb 10, 2022, at 1:32 AM, Yasmin97 <yasmin...@gmail.com> wrote:

Hi Wendy, 

That is a pretty huge difference in marginal likelihoods. If I were you, I would double check my data sets or the analysis settings? I have never seen such a massive difference between models! At most, I have seen a difference of perhaps a few hundred between two models (however, I have never compared a strict to a random clock). 

Do you have the output from the generalised-stepping-stone analysis? 

Best,

Yasmin 

On Thursday, February 10, 2022 at 7:19:58 AM UTC+11 weedy23 wrote:

-- 
You received this message because you are subscribed to the Google Groups "beast-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to beast-users...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/beast-users/9c5bce63-7216-4e70-890a-a55f5fe4867en%40googlegroups.com.

[David Swofford]

I guess I probably should have written "a *Metropolis-Hastings* MCMC sampler *that has not been designed to deal cleanly with out-of-bounds proposals* may go haywire. Technically, it should still converge "eventually" but low acceptance probabilities may cause the chain to become stuck for a very, very long time (possibly longer than the machine can be kept running).

Dave

To view this discussion on the web visit https://groups.google.com/d/msgid/beast-users/EA4E600A-DDA5-46CE-88B3-3857A72F6C72%40duke.edu.

[weedy23]

Hi Dave,

Thank you for your reply. I guessed as much, however I am still at a loss for what could have gone wrong. This is the same analysis I have performed before (without model selection applied) and so there doesn't seem to be a problem with my input data. Potentially my priors, however they seem pretty reasonable and the Tracer outputs from these analyses doesn't seem to identify any major issues. Is there any more information I can post to help identify the issue?

Many thanks for your help,

Sarah

[weedy23]

Hi Yasmin,

Thanks for your help. Which output are you talking about sorry? Something different to the output I posted above?

[Yasmin Asar]

Hi Sara,

Yes, I see you're using BEAST (v1.10.4?) I use BEAST2 and the output looks a little different, sorry for not picking that up. 

See this paper: Improving the Accuracy of Demographic and Molecular Clock Model Comparison While Accommodating Phylogenetic Uncertainty (Baele et al. 2011)

" To ameliorate these issues, we simplified the GTR model to an HKY (Hasegawa, Kishino, and Yano) model, fixed the base frequencies to the empirical base frequencies and most importantly replaced the improper uniform prior on the mean rate in the UCLD model with a diffuse gamma prior...This resolved the apparent mixing issues, yielding proper posterior and prior distributions and consistent model ordering according to PS and SS (table 2). It therefore remains a crucial part of any MCMC analysis to check the MCMC chain for adequate mixing and provide proper priors for all the model parameters if one wishes to estimate marginal likelihoods. "

Perhaps it would be useful for you to set an upper limit on the clock rate prior (if currently, it is unbounded to infinity)? I had a similar issue in the past, but we fixed it by making sure the 'ucldmean' was a uniform distribution fixed bounded between 0 and 1, with an initial value of 0.01. Perhaps you could double check the value for your clock rate prior? 

Also see Duchene et al. (2020): "We specified proper priors on all parameters, which is essential for accurate estimation of log marginal likelihoods (Baele et al. 2013). In our heterochronous analyses the prior on the evolutionary rate had a uniform distribution bounded between 0 and 1."

Best,

Yasmin 

[weedy23]

Thanks Yasmin, that is very helpful. Yes I'm using BEAST v.1.10.4. Currently my clockrate is using the default CTMC rate reference, which I had read somewhere was usually appropriate, but yes maybe not for MLE. I will try specifying a proper prior on this and see if I see any improvement. 

Thanks again for your help, much appreciated.

[Yasmin Asar]

No problem at all Sara! Let us know how it goes!

Best,

Yasmin

[Sarah Wells]

Hi Yasmin, and anyone else!

So I changed the clock.rate prior to uniform as you suggested but it didn't change anything unfortunately. I have been fiddling around with various other things but I am still getting huge (albeit slightly smaller) log likelihoods of -4e70 when using the random local clock.

I wonder if it could be something to do with a poorly formed XML file when generalised stepping stone estimation is selected, because I have been having issues with BEAUTI generated files (in some versions there were incomplete MLE blocks generated).

For example, I have been getting errors and terminations during the MLE process stating

"Error parsing '<compoundLikelihood>' element with id, 'null': An element (yuleReference) which is not a likelihood has been added to a compoundLikelihood element".

This seems to refer to the MLE settings block where it keeps putting in <yuleModel idref="yuleReference"/> in the working prior where I think there is supposed to be the speciationlikelihood? In any case, I changed this line to <speciationLikelihood idref="speciationReference"/> (see line 1035 in attached file) and MLE then ran (but with the huge log likelihood end result), but I'm now wondering if this was correct.

I am also getting a warning stating

"WARNING: Likelihood component, null, created but not used in the MCMC".

These XMLs are generated by BEAUTI and unmodified by me, so I'm not sure whether its a bug or something wrong in my settings. Would someone please be able to have a look at the XML attached (with my above change made on line 1035) to see if they can spot what is going wrong and if anything is out of place? 

I’m stumped.

--
You received this message because you are subscribed to a topic in the Google Groups "beast-users" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/beast-users/ZjUM3wySs9o/unsubscribe.
To unsubscribe from this group and all its topics, send an email to beast-users...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/beast-users/95dfcfe9-409a-4a4f-8a2d-1113d60435cdn%40googlegroups.com.

[Yasmin Asar]

Hi Sara, 

I notice that in your .XML file, you have a chain length of 1,000,000. But, in your stepping-stone sampling run you have 50 steps with a chain length of 100,000, totalling a chain length of 5,000,000. Just note that the length of the GSS run should be the same as the total MCMC chain length. To do this, you would reduce the number of steps to 10, since 10 x 100,000 = 1,000,000. I am really not sure that this is causing your issue, however, since you get such crazy log likelihoods. I'm definitely quite stumped with this! 

Would you consider using BEAST2 to run your data? I have not had any issues with the MODEL SELECTION package that is provided, and have quite smoothly run a lot of generalised stepping-stone analyses there! Would recommend if you're getting buggy issues with BEAST1. 

I hope you manage to get through this!

Yasmin. 

[weedy23]

Hi Yasmin,

This XML file was a test, because it was taking so long to run and I kept getting silly likelihood values, I reduced the length of the MCMC and GSS when I was trying to get it to work. The original runs actually had an MCMC of 100M, and 100 steps of 1M length in the GSS. So yeah it is not that causing the huge log likelihoods. 

Thanks for suggesting BEAST2, I may switch if I really can't work it out. Does BEAST2 have GSS as well or is it better to use a different method?

Many thanks,

Sarah

[Yasmin Asar]

Hi Sara,

Yes, BEAST2 has GSS! You can definitely use it. You could try nested sampling, but it is very tricky to get the parameters right with nested sampling.

Best,

Yasmin