Wondering about model comparison, PS/SS in BEAST v1.71.

[Giap Nguyen]

Dear BEAST users!

 

Recentlty, I run analysis about model comparison according to Guy Baele et al. (MBE, 2012). I also obtained .xml file as supplement data and tried my self on BEAST v.1.7.1. However, I do not have a general idea when the MCMC will terminate (now it is said that "attemping theta = 4.19... e-4).

Would you please show me a little more detail on model comparison by BEAST?

 

Thank you very much!

 

Giap Nguyen

[Guy Baele]

Dear Giap,

Path sampling and stepping-stone sampling run a series of power posteriors, starting at 1.0 and ending at 0.0.
Hence, the closer you get to 0.0, the closer your analysis is to terminating.

The path from 1.0 to 0.0 is not uniform however (although you could do this), but follows a beta distribution.
This means that half of the theta values will be below 0.1 and that more theta values are traversed near 0.0.
Any way, given the value you mention, I'm sure the analysis will finish soon enough.

Best regards,
Guy

Op maandag 21 mei 2012 09:17:24 UTC+2 schreef nvgiap het volgende:

[Giap Nguyen]

Dear Dr. Guy Baele

I would like to thank you very much for your reply.

Sincerely

Giap Nguyen

[Giap Nguyen]

Dear Dr. Guy Baele

 

I would like to ask about several things that I concern about model selection using PS/SS.

 

1. Do I need to seperate the code for "marginalLikelihoodEstimator" and "pathSamplingAnalysis/steppingStoneSamplingAnalysis" in the two different .xml file to run with BEAST?

 

2. What parameter do I need to use (in the following example) to compare between models?

 

Thank you very much!

Sincerely yours

Giap Nguyen

 

1 file(s) found with marginal likelihood samples
PathParameter MeanPathLikelihood
0.95087 -6452.1
0.95409 -6456.1
0.95732 -6457.6
0.96056 -6456.3
0.96380 -6454.1
0.96705 -6455.8
0.97031 -6452.6
0.97358 -6457.6
0.97686 -6451.6
0.98014 -6450.3
0.98343 -6454.1
0.98673 -6445.4
0.99003 -6450.3
0.99335 -6446.6
0.99667 -6452.8
1.0000 -6450.4

log marginal likelihood (using path sampling) from pathLikelihood.delta = -317.03130576335496

Inner area for path parameter in (0.95409,0.99667) = -296.24

 

1 file(s) found with marginal likelihood samples
PathParameter MaxPathLikelihood
0.95087 -6435.4
0.95409 -6436.4
0.95732 -6437.1
0.96056 -6436.3
0.96380 -6432.9
0.96705 -6439.8
0.97031 -6429.0
0.97358 -6436.1
0.97686 -6429.5
0.98014 -6430.4
0.98343 -6426.5
0.98673 -6429.2
0.99003 -6429.3
0.99335 -6425.0
0.99667 -6429.1
1.0000 -6432.7

log marginal likelihood (using stepping stone sampling) from pathLikelihood.delta = -317.0283600379283 

[Guy Baele]

Dear Giap,

1. I'm not sure what you mean by 'separating the code'. The marginalLikelihoodEstimator element runs a series of power posteriors between posterior and prior and, in doing so, collects samples from each of these power posteriors. What it does not do however, is provide an estimate of the log marginal likelihood. This is done by the pathSamplingAnalysis and/or steppingStoneSamplingAnalysis elements, which take the collection of sampels, and estimate the log marginal likelihood from these samples.

2. You need to use the 'final result', i.e. the estimate of the log marginal likelihood, to compare models. For path sampling, this is:

log marginal likelihood (using path sampling) from pathLikelihood.delta = -317.03130576335496

And for stepping stone sampling, this is:

log marginal likelihood (using stepping stone sampling) from pathLikelihood.delta = -317.0283600379283 

The other output that you see just shows part of the calculations that are being done, but you can't use it to perform model selection.

Cheers,

Guy

Op vrijdag 8 juni 2012 04:01:10 UTC+2 schreef nvgiap het volgende:

Dear Dr. Guy Baele

 

I would like to ask about several things that I concern about model selection using PS/SS.

 

1. Do I need to seperate the code for "marginalLikelihoodEstimator" and "pathSamplingAnalysis/steppingStoneSamplingAnalysis" in the two different .xml file to run with BEAST?

 

2. What parameter do I need to use (in the following example) to compare between models?

 

For path sampling

[Giap Nguyen]

Thank you very much for your reply!

 

Sincerely

Giap Nguyen

2012/6/8 Guy Baele <bael...@gmail.com>

--
You received this message because you are subscribed to the Google Groups "beast-users" group.

To view this discussion on the web visit https://groups.google.com/d/msg/beast-users/-/Vz-dh72sFWsJ.

To post to this group, send email to beast...@googlegroups.com.
To unsubscribe from this group, send email to beast-users...@googlegroups.com.
For more options, visit this group at http://groups.google.com/group/beast-users?hl=en.

[Roberto Gomez]

Dear Guy,

Thanks for all the help.  I think Giap may have been confused by the same issue that has me confused.  I think he was asking about was whether the other commands need their own xml files separate from the input file that specifies SS/PS run details. 

I'm performed several BEAST analyses, and I'm interested in using PS/SS to do hypothesis testing in a Bayesian framework since I want to do the method that the community has agreed upon as being the currently best practice.  I had no trouble running the BEAST analyses, and all the ESS values are high.  Where I get confused are the specifics of calculating the Bayes factor from the marginal likelihood to have some statistical testing to see whether one hypothesis (clade A is x old) fits the data better than another hypothesis (clade A is y old) by .....   

I recently learned that through BEAUti one can include a MLE mcmc run in the input file.  I selected 100 steps and 15 million generations, and I left everything else the same as my previous BEAST analyses. 

What I am still uncertain of is how to use the output values from the PS/SS to compare models.  Is it simply the marginal likelihood of model 2 subtracted from the marginal likelihood from model 1, and at that point is it actually Bayes factor?  Is the calculation a ratio between models?  Also, how do you decide (even if it's arbitrary) if the MLE mcmc has run long enough? 

I've been really stumped by this, and it's mostly that I enjoy phylogenetics and want to make sure I'm doing the right thing especially as a student.

Cheers,
Roberto 

[Guy Baele]

Dear Roberto,

The Bayes Factor is a ratio of two marginal likelihoods and hence the log Bayes Factor is the difference between two log marginal likelihoods.
BEAST estimates the log marginal likelihood (1 per BEAST run), so yes, you have to subtract the log marginal likelihood for model A from the log marginal likelihood for model B.

Your settings may have been quite a bit exaggerated if I understand them correctly, i.e. you're running 100 power posteriors and each power posterior runs for 15 million iterations? That's a total of 1.5 billion iterations (more if you count the burn-in that's not being logged to file), which might take months to compute.

Trying 100 power posteriors, with 1 million iterations per power posterior would be a good start in my opinion, although it depends on the size of your data set and the kind of data you're using.

To know if it has run long enough, there's not really a statistic you can compute. It's often suggested to run the MLE calculations again with different settings. For example, if it's feasible, try running it twice as long (by for example doubling the number of power posteriors). If the marginal likelihood estimates are very close to the ones you obtained in the first run, you're good and you don't have to look any further.

Best regards,

Guy

Op dinsdag 14 oktober 2014 02:27:17 UTC+2 schreef Roberto Gomez:

[Mariana López]

Deary Guy,

a few years later, I have two questions:

I used Beast2 to calculate SS/PS to compare clock models. Reading the previous messages, I have two questions:

1. I used beast2, MODEL_SELECTION package, and I get the result I show below. I don't understand if this metho is a combination of PS and SS because I don't see any different result for SS or PS as mentioned before. In the tutorial I used (https://github.com/BEAST2-Dev/beast-docs/releases/download/v1.0/BFD-tutorial-2018.zip) there is no differentiation between SS and PS.. can you clarify this, pelase?

2. Based on your previous comments I understand that my results are ok, I mean theta reach 0 and half of the steps are below 0.1 (I run 100 steps and 1 million iterations). Am I right?

Thanks in advance

Mariana

Step        theta         likelihood   contribution ESS
0            1            -52778.4083  0            5.0409      
1            0.9667       -52754.2809  -1754.0963   4.3971      
2            0.9342       -52784.3059  -1713.8696   5.2674      
3            0.9025       -52813.6066  -1674.1577   5.8041      
4            0.8716       -52808.1096  -1633.2276   4.1962      
5            0.8413       -52863.2015  -1594.72     3.8121      
6            0.8119       -52825.314   -1555.0452   4.9211      
7            0.7831       -52851.6474  -1517.2074   5.0824      
8            0.7551       -52855.1063  -1479.5339   4.4086      
9            0.7278       -52859.741   -1442.1288   4.5656      
10           0.7012       -52911.666   -1406.6571   5.2055      
11           0.6753       -52868.061   -1369.2344   4.2684      
12           0.6501       -52973.1052  -1335.7664   4.7848      
13           0.6255       -52966.2085  -1299.7253   3.8544      
14           0.6016       -53004.6886  -1265.707    3.942       
15           0.5783       -52973.1948  -1232.0777   4.6998      
16           0.5557       -53000.364   -1198.8901   9.0795      
17           0.5337       -53013.7316  -1165.6147   7.3451      
18           0.5123       -53038.1137  -1133.3857   4.5551      
19           0.4915       -53078.3333  -1101.9777   4.749       
20           0.4713       -53110.7391  -1071.3687   9.9913      
21           0.4517       -53181.0345  -1040.953    3.9717      
22           0.4327       -53197.1496  -1010.6228   4.0507      
23           0.4142       -53251.7683  -981.3079    3.4991      
24           0.3964       -53307.0043  -952.6848    4.1166      
25           0.379        -53236.9157  -922.882     6.101       
26           0.3622       -53303.0488  -895.3052    5.787       
27           0.3459       -53427.1047  -868.7834    4.7305      
28           0.3302       -53420.4209  -841.1066    4.3036      
29           0.3149       -53461.4777  -814.6603    10.2296     
30           0.3002       -53451.3325  -788.0585    16.852      
31           0.2859       -53507.7997  -762.5578    18.4916     
32           0.2721       -53554.3459  -737.6001    22.6067     
33           0.2588       -53652.2539  -713.528     14.7753     
34           0.246        -53769.6357  -690.0637    6.855       
35           0.2336       -53752.5891  -665.6962    9.9359      
36           0.2217       -53806.4237  -642.2418    6.5249      
37           0.2101       -53981.5895  -620.8488    5.0816      
38           0.1991       -54119.4591  -598.8003    9.989       
39           0.1884       -54151.7115  -577.453     16.7516     
40           0.1781       -54214.863   -555.9716    11.6391     
41           0.1683       -54395.3736  -536.1378    16.2526     
42           0.1588       -54483.1218  -515.7671    24.1206     
43           0.1497       -54624.829   -495.8211    5.0096      
44           0.141        -54783.3041  -477.5059    19.0801     
45           0.1326       -55083.0566  -459.6682    8.6631      
46           0.1246       -55262.1959  -442.0888    34.3322     
47           0.1169       -55588.5148  -424.7315    9.0949      
48           0.1096       -55906.4812  -407.2495    9.7335      
49           0.1026       -56152.0884  -392.0756    4.9588      
50           0.0959       -56647.957   -376.734     6.9896      
51           0.0895       -56905.4168  -360.6965    13.7017     
52           0.0835       -57207.0476  -346.3766    11.5282     
53           0.0777       -57642.1345  -330.9277    8.5693      
54           0.0722       -57962.9446  -315.1279    11.0017     
55           0.067        -58359.3742  -302.0632    13.3012     
56           0.0621       -58589.4609  -288.9549    39.1875     
57           0.0574       -59055.9971  -275.9583    31.0294     
58           0.0529       -59229.5821  -261.4976    27.3022     
59           0.0488       -59527.3754  -248.6104    50.5859     
60           0.0448       -59867.2715  -235.8092    64.4773     
61           0.0411       -59976.203   -222.5523    49.4464     
62           0.0376       -60231.8624  -210.2166    55.1265     
63           0.0343       -60388.2298  -197.9402    51.9612     
64           0.0312       -60511.8801  -185.9035    46.1675     
65           0.0284       -60681.912   -174.4014    60.8029     
66           0.0257       -60853.4665  -163.3396    16.236      
67           0.0232       -60932.204   -152.3811    93.5128     
68           0.0208       -61134.1305  -142.1413    52.1024     
69           0.0187       -61239.0712  -132.0397    24.9337     
70           0.0167       -61252.9821  -122.2026    19.7047     
71           0.0149       -61492.307   -113.1409    8.0914      
72           0.0132       -61459.2749  -104.1134    79.6281     
73           0.0116       -61649.6163  -95.783      32.0085     
74           0.0102       -61756.2476  -87.695      15.3042     
75           0.0089       -61885.9112  -80.0649     25.2909     
76           0.0077       -61927.1535  -72.699      26.1291     
77           0.0066       -62120.8048  -65.897      75.4385     
78           0.0057       -62384.9946  -59.4824     10.0911     
79           0.0048       -62258.9542  -53.1463     55.943      
80           0.0041       -62499.1029  -47.4614     81.6126     
81           0.0034       -62446.1115  -41.9453     57.8794     
82           0.0028       -62946.2565  -37.1352     22.1408     
83           0.0023       -63568.1711  -32.6654     10.7074     
84           0.0019       -63354.2424  -28.1587     48.687      
85           0.0015       -63788.7097  -24.2507     64.1554     
86           0.0012       -63956.7125  -20.5686     15.4958     
87           0.0009       -64211.4242  -17.2334     61.1579     
88           0.0007       -64849.6036  -14.3298     55.9232     
89           0.0005       -65917.8856  -11.7681     59.9161     
90           0.0003       -67000.2793  -9.463       112.933     
91           0.0002       -67726.2846  -7.3658      69.3468     
92           0.0001       -69703.8164  -5.6592      100.608     
93           0.0001       -71645.3354  -4.1443      60.207      
94           0            -74026.1733  -2.8922      67.5556     
95           0            -75174.1281  -1.8565      73.8221     
96           0            -77715.7186  -1.0715      70.3122     
97           0            -77941.9392  -0.4992      38.3229     
98           0            -79138.5579  -0.1597      24.496      
99           0            -79712.8712  -0.0178      87.5777     
sum(ESS) = 2685.9576  

marginal L estimate = -53863.1049830501