Bayes factors for comparing variances

This post is adapted from my Observable notebook. The package results discussed below are those recorded in that notebook.

Introduction

Dablander et al. (2024) outline a method for computing Bayes factors for variance comparisons. Section 5.2 gives an example of comparing a single sample variance with a reference variance, using both an undirected and a directed hypothesis.

The example is as follows:

The EPA is testing a new device for measuring PCB concentrations in fish. The measurements of the device are required to have a variance of less than 0.1 (). Therefore, . Seven readings are taken from a sample of fish, and the sample standard deviation is . This gives a sample effect size

Dablander et al. compare two hypotheses:

The Bayes factor can be calculated with the accompanying R package as follows:

> bfvartest::onesd_test(n=7, s=.22, popsd=sqrt(0.1), u = 0.5, logarithm=FALSE, alternative_interval=c(1,Inf))
[1] 0.509796

This matches the value reported in the paper.

Dablander et al. also reported the undirected Bayes factor. That is, the Bayes factor from comparing the following hypotheses:

This can be calculated with the accompanying R package as follows:

> bfvartest::onesd_test(n=7, s=.22, popsd=sqrt(0.1), u = 0.5, logarithm=FALSE)
[1] 0.4104968

Again, this matches the value reported in the paper.

Bayesplay for variance Bayes factors

The Bayesplay package also has the functionality to compute Bayes factors for differences between variances. These are parametrised slightly differently in Bayesplay compared with bfvartest, but they work in essentially the same way. In Bayesplay, the likelihood is parametrised in terms of u, which is the variance ratio of the sample variance against the reference variance (), and a df of n-1.

Inputting the values from Example 5.2, we get the following:

Because the likelihood is in terms of variance, and not precision, the bounds of the truncated prior are also specified in the opposite direction.

With these differences in mind, we can specify the models in the Bayesplay modelling language as follows (both directed hypotheses are included, because the results will become useful below):

model = `
    analysis R {
      likelihood variance_ratio(u = 0.484, df = 6)

      hypotheses {
        0: point(1)
        1: beta_prime(0.5)
        -: beta_prime(0.5) | R < 1
        +: beta_prime(0.5) | R > 1
      }
    }
`

Again, note that because the likelihood is in terms of variance rather than precision, above is the same as the reported by Dablander et al.

If we look at the results produced by Bayesplay, we can see that the undirected Bayes factor is the same as that produced by bfvartest and the value reported in the paper: .

However, if we look at the directed Bayes factor, the result produced by Bayesplay is different. For the directed Bayes factor, Bayesplay gives the result: .

Which result is correct?

First, we can look at the two directed hypotheses and whether they’re consistent with the undirected hypothesis.

We should expect that .

Checking the numbers produced by bfvartest we get:

We can also check the numbers produced by Bayesplay:

We can get a clue as to what is going wrong with bfvartest if we change the model we use with Bayesplay slightly. Instead of defining as (equivalent to bfvartest with ), we define as follows:

This gives the following model.

model2 = `
    analysis R {
      likelihood variance_ratio(u = 0.484, df = 6)

      hypotheses {
        0: point(1)
        1: beta_prime(0.5)
        -: beta_prime(0.5) | R < 10
      }
    }
`

With this new model definition, we obtain a , which accords with the value reported in the paper.

We can also update the prior, to use the informed prior reported in the paper. In the Bayesplay modelling language, this gives the following model:

model3 = `
    analysis R {
      likelihood variance_ratio(u = 0.484, df = 6)

      hypotheses {
        0: point(1)
        1: beta_prime(2.16)
        -: beta_prime(2.16) | R < 10
      }
    }
`

With this model definition, we get a result of , which again accords with the result in the paper.

These findings taken together suggest that there is an error with bfvartest where the directed priors are truncated relative to rather than being truncated at 1.

To confirm this, we can modify the model again. We can change our data slightly so that we’re measuring our data on a slightly different scale, for example, parts per 100,000 rather than parts per 10 million. For Bayesplay the model is unchanged, so the results are unchanged.

For bfvartest we’d update the function call as follows:

> bfvartest::onesd_test(n=7, s=2.2, popsd=3.2, u = 0.5, logarithm=FALSE, alternative_interval=c(1,Inf))
[1] 7.897717e-05

We should not expect our Bayes factor to change just because we change our units of measurement. This again points to an error with bfvartest.

We can attempt to replicate this result with Bayesplay by truncating the prior at rather than 1. Doing this produces a of 7.6e-5, which is the same order of magnitude as the result of bfvartest.

The issue with bfvartest

If we have a look at the code for bfvartest, specifically the function .compute_logml_restr_k1 in helpers.R, we can see where the error is. If we correct this function, then bfvartest produces the correct results.

Try Bayesplay yourself

Below you can try specifying a model using a web interface to Bayesplay.

The default values are the same values as the example in Section 5.2 of Dablander et al. (2024). You can change the sample standard deviation, the reference standard deviation, the sample size, and the prior alpha. You can also set where to truncate the prior for each of the two directional hypotheses.

The Bayes factor is shown for the two-sided comparison and each of the two directional hypotheses.

References

Original Observable notebook, first published on 23 August 2026. Calculations use Bayesplay JS at revision 0ffabea3.