Ryan Schott sends in this question, which he has kindly agreed to be posted here:

"I, along with other members of Belinda Chang's lab at the University
of Toronto, have been using FUBAR, along with other HYPHY and PAML
methods, to study variation in dN/dS across sites. In some datasets,
but not all, we are finding that FUBAR has elevated dN/dS ratios
relative to both REL and M8/M2. This is most noticeable in an elevated
baseline when the values are plotted on a graph. An example of this
can be seen in the attached figure. It seems odd to us the FUBAR and
REL should differ so much, but REL and M8 be similar. With a dataset
we are currently working on the difference between FUBAR and REL/M8 is
even larger and we are currently at a loss as to why that is the case.
With other those the values are similar between all models.

I am hoping you might have some ideas on why FUBAR and REL might
differ substantially in (baseline) dN/dS values. Thanks for your time!"

My response:

The first thing to bear in mind is that none of these methods (all of these are random effects methods, but this is also true for fixed effects methods) are designed or expected to always give reliable point estimates of dN/dS. What has been evaluated for all of these methods is the statistical performance of the hypothesis tests, e.g. is dN/dS>1? When we see a reported dN/dS point estimate of, say, 4.2, this does not give us a confidence interval (or other measure of how certain/uncertain the point estimate is). Sometimes (when the data set is large and informative) we may be pretty sure that dN/dS really is close to 4.2, but sometimes (when the data set is small or uninformative) we might get the same estimate when all we really know is that dN/dS is somewhere between 0.5 and 10. So instead of relying on point estimates, these methods use posterior probabilities to decide whether the evidence strongly favours dN/dS>1.

Next, consider the main difference between FUBAR and the other REL methods: FUBAR uses 400 rate categories, whereas the other methods use only a handful. The effect of this is that FUBAR has much more flexibility to obtain estimates anywhere between say 0 and 10, whereas the other methods are more likely to get estimates close to the estimated category values (the estimates are weighted sums of all the categories so they could be in-between, but often one category dominates and then the estimate is close to the dN/dS value of that category). In your graph, we see this in two places: 1) For sites under purifying selection, REL and M8 have many estimates close to zero, presumably because REL placed one of its three categories there and M8 placed most of the weight of its beta distribution there. FUBAR, on the other hand gets intermediate values (closer to 1) for many of the sites, because it also has categories at those intermediate points. 2) For sites under positive selection, both REL and M8 estimate the dN/dS values at most of the positively selected sites to be roughly equal to each other - the estimates are dominated by where the single positive selection category was placed. FUBAR estimates different dN/dS values at different positively selected sites - some larger and some smaller than the REL/M8 values. This behaviour of the non-FUBAR methods is an undesirable artifact of having a small number of rate categories - we have no reason to believe that dN/dS should be similar at different positively selected sites.

Finally, consider the sites that are actually marked as significant (high posterior probability of being under positive selection). We see several examples identified by REL and/or M8 but not by FUBAR. Here, the REL/M8 dN/dS estimates have been dragged up toward the position of the positive selection category, which could cause false positives - this is an example of the effect we demonstrated in Fig 3 of the FUBAR paper. Interestingly, in this example we see no sites identified by FUBAR but not by REL/M8 (this need not have been the case) - which underscores that FUBAR is not spuriously generating inflated dN/dS estimates.

Thanks Konrad for your very detailed answer.

The difference in the number of rate categories was my best guess (after ruling out variation in dS), but didn't really make sense to me since the point estimates did not seem constrained to particular rate categories. Your explanation of how the point are calculated as weighted sums of all the categories really clears this up.

What is interesting to note is in the datasets where we see this elevated baseline in FUBAR compared to REL/M8 the distribution lanscape of dN and dS is very complex (much more variation in rates than could be accounted for by REL or M8), but in thos that are similar between FUBAR and REL/M8 the rate distribution is much less complex. This really highlights the power of FUBAR to analyze complex datasets.