Note: Student made a few mistakes when entering the data from Cushny and Peeble study (documented in (2)). He mislabeled the columns and omitted the observation for Patient 11.

Student (9) analyzed the data from a Bayesian perspective and considered only the two hypotheses H_l and H_m (i.e. a one-sided test). He computes probabilities consistent with using the prior

Note: The prior Student uses is improper, which perhaps requires some discussion. For the special case of one-sided testing, improper priors like this can be justified though the prior assumes conditions that are in a sense least favorable to the null hypothesis (see (11) and (12, p. 147–148)).

Here’s some python code that reproduces the values from Student’s paper:

import numpy as np
# Note: Student only considers the first 10 entries from the data set
deltas = np.array((
     -1.2, -2.4, -1.3, -1.3,  0. , -1. , -1.8, -0.8, -4.6, -1.4))
n = len(deltas)

t = np.mean(deltas) / (np.std(deltas, ddof=1) / np.sqrt(n))


# Compute the posterior probability for the hypothesis that
# L-hyoscyamine HBr is less effective than L-hyoscine HB
import scipy
P_H_alt = scipy.stats.t(n-1).cdf(-t)
print('Probability L-lyoscine HB is more effective:', P_H_alt)

Probability L-hyoscine HB is more effective: 0.9985835549013079

Rounding will give us the same result that Student reports.

To make our results more comparable to Student’s, we’ll adjust the values for the n=10 entry in Table 2 to exclude the equals hypothesis, giving us

or an odds of about 119 to 1. As expected, EEIBFs give a smaller probability than Student since the method is design to provide a more reasonable prior probability for the null hypothesis than what would correspond to a lower bound.

The building blocks for objective Bayesian hypothesis testing are the Bayes factor and objective priors.

Let M_1 and M_2 denote two competing models with parameters θ_1 and θ_2 and proper priors π_1 and π_2. Suppose we then observe data x. Then the posterior probability for one of the models is given by

The Bayes factor, B_ji, is defined as the ratio

If the models are given equal prior probability, then B_ji is the ratio of the posterior probabilities for the two models.

Suppose we observe 5 successes and 4 failures from a binomial distribution with unknown probability of success parameter p. Now imagine that we consider two models corresponding to the hypotheses that p1/2. If we use Jeffreys prior, then the priors for the two models after renormalization are

and the Bayes factor is given by

Let’s now consider the normal distribution. The reference prior for the normal distribution with the mean, µ, as the parameter of interest and the variance, σ², as the nuisance parameter is given by

(see example 10.5 of (17)). In this case, π is improper so we can’t use it directly to compute Bayes factors. Improper priors are unnormalized and only specified proportionally. Any multiple c x π(µ, σ²) is also a valid version of the improper prior; and, hence, computing Bayes factors with improper priors gives us a result that’s dependent on the ratio of arbitrary constants.

One way to work around this issue, the training sample method (10), is to use a minimal subset of the data to form proper priors and then the remaining data to compute the Bayes factors.

Note: By minimal, I mean the smallest subset of observations such that the posterior using the prior π is proper and well defined.

Suppose we observe the values

from a normal distribution with known variance σ²=1 and we want to test the two hypotheses, µ 0. For the objective priors

we need a single observation to form proper distributions. Using the first observation, gives use the new priors

where φ and Φ represents the PDF and CDF of the standard normal distribution. We can then compute the Bayes factor

Clearly, the Bayes factor computed in Example 2 will depend on which observation is chosen to “train” the prior. As an improvement, we can average the result over all possible minimal training subsets, and this leads us to what Berger and Mortera (8) call encompassing arithmetic intrinsic Bayes factors (EIBFs).

Suppose H_1 and H_2 are two hypotheses to be tested and let H_0 denote a hypothesis that encompasses both H_1 and H_2.

For example, if H_1 and H_2 denote hypotheses that the mean of a normal distribution is less than zero and greater than zero, then a naturally encompassing hypothesis is that the mean is unrestricted.

The encompassing arithmetic intrinsic Bayes factor is defined as

where B_ji^N denotes the Bayes factor computed using the improper prior and x(l) denotes a minimal subset of the observed data and l enumerates over all possible minimal subsets.

EIBFs can work well with enough observations; but for smaller data sets, the term

can be numerically unstable. To make the factors more general purpose, Berger and Mortera recommend replacing the sum with expected values which gives us the expected encompassing intrinsic Bayes factor (EEIBF),

where X(l) denotes a minimal vector of independent identically distributed random variables large enough to make the posterior proper and distributed according to the model for the encompassing hypothesis H_0 with parameters that maximize the likelihood of the observed data x.

Suppose we observe the same values from Example 2,

and assume the data was sampled from a normal distribution with unknown variance and we wish to test the two hypotheses, H_1: µ 0. For the objective priors,

we need two observations to form a proper posterior distribution. The encompassing hypothesis, H_0, is that the mean is unrestricted and the maximum likelihood parameters are given by

Thus, the expected encompassing intrinsic Bayes factor is given by

where X_1 and X_2 are independent normally distributed variables with mean µ_ML and variance σ²_ML.

We can use the package bbai to do the computation.

import numpy as np
y = np.array((1.89, 0.52, 1.10, 2.36, 1.99, -0.85, 1.07))


from bbai.stat import NormalMeanHypothesis
result = NormalMeanHypothesis(with_equal=False).test(y)
print('B_12^EEI =', result.factor_left / result.factor_right)

B_12^EEI = 0.07224015471356342

Berger and Mortera write:

I’ll add that for the case of normal mean testing with unknown variance, there’s no need to resort to simulation. With appropriate quadrature rules and interpolation in Chebyshev nodes, we can accurately compute EEIBFs deterministically and efficiently. Perhaps I’ll go more into the details in another blog post; but for now, there’s this Jupyter notebook that shows how to do a step-by-step validation of an efficient deterministic implementation: https://github.com/rnburn/bbai/blob/master/example/18-hypothesis-eeibf-validation.ipynb.

Q1: What other options are there for objective Bayesian hypothesis testing?

As mentioned in the introduction, Sellke and Berger developed equations that allow you to derive a lower bound for the posterior probability of the null hypothesis, and Clare (16) makes improvements on Sellke and Berger’s work to give tighter bounds.

As another approach, Berger and Mortera developed intrinsic priors that asymptotically give the same answer as various default Bayes factor methods such as EEIBF, and they suggest that these intrinsic priors might be used in place of default Bayes factors:

Q2: But if there are multiple reasonable methods, doesn’t that make the EEIBF approach subjective, whereas P-values are an objective approach to hypothesis testing?

If the goal is to derive an objective answer to a hypothesis test, then P-values are a red herring.

A P-value is determined not just by the model and observed data, but also by experimental intent (22). The same observation and model can result in remarkably different P-values. Berger and Berry (13) provide this example to demonstrate.

Suppose a research studies a coin for bias. The researcher flips the coin 17 times. Heads comes up 13 times and tails comes up 4 times. Let θ represents the probability of heads and suppose the researcher is doing a standard P-value test with the null hypothesis that the coin is not bias, θ=0.5. What P-value would they get?

If the researcher intended to flip the coin 17 times, then the probability of seeing a value less extreme than 13 heads under the null hypothesis is given by summing binomial distribution terms representing the probabilities of getting 5 to 12 heads,

which gives us a P-value of 1–0.951=0.049.

If, however, the researcher intended to continue flipping until they got at least 4 heads and 4 tails, then the probability of seeing a value less extreme than 17 total flips under the null hypothesis is given by summing negative binomial distribution terms representing the probabilities of getting 8 to 16 total flips,

which gives us a P-value of 1–0.979=0.021.

One can argue whether any statistical method can be truly objective or to what extent objectivity should even be the goal; but methods based on objective Bayesian inference such as EEIBF have at least as strong claim to objectivity as P-values.

Q3: Does rejection of P-value hypothesis testing make one a Bayesian?

We shouldn’t equate frequentism with P-value hypothesis testing, and many of the results from objective Bayesian inference can be adapted to use frequentist terminology.

In (18), for example, Berger considers the thought experiment of whether Fisher, Jeffreys, and Neyman’s views on hypothesis testing could be reconciled and he shows how Bayes factor values could be interpreted by a frequentist as likelihood ratios and used to produce conditional error probabilities.

When Fisher initially advocated for P-value testing, the method was met with skepticism and produce a lively debate. Harold Jeffreys famously criticized the use of P-values:

We can certainly appreciate the problem Fisher was trying to solve with P-values. At the time, common practice was to use Bayesian statistics (i.e. inverse probability) with a default uniform prior. Such an approach will, of course, produce results that depend arbitrarily on the parameterization of the model (20, 21). But as Berger and Berry show (13), P-values can be just as arbitrary, and the frequent misinterpretation of P-values as probabilities shows us that statistics users naturally want to interpret results as a measurement of degree of belief.

Default Bayes factors such as EEIBF fix the problems of early Bayesian practice by replacing naive use of the uniform prior with reference priors while keeping the goal of delivering results in terms of degree of belief. While no such method can claim to be objective in the purest sense (I doubt a completely objective approach to hypothesis testing exists), they can make a reasonable claim to objectivity and they can produce good results under very general circumstances.

(1): Cushny, A. R. and Peebles, A. R. The action of optical isomers. 11: Hyoscines’, Journal of Physiology, 32, 501–510 (1905). Since the paper is from more than 95 years ago, the data set can be consider to be in the public domain.

(2): Senn S, Richardson W. The first t-test. Stat Med. 1994 Apr 30;13(8):785–803. doi: 10.1002/sim.4780130802. PMID: 8047737.

(4): Bernstein, P. (1998). Against the Gods: The Remarkable Story of Risk. Wiley.

(12): Berger, J. (1985). Statistical Decision Theory and Bayesian Analysis. Springer.

(14): “View of Michigan Asylum for the Insane, Kalamazoo” Map of Kalamazoo Co., Michigan. Philadelphia: Geil & Harley, 1861. Library of Congress.

(15): Fisher, R. A (1925). Statistical Methods for Research Workers.

(20): Fisher, R. (1930). Inverse probability. Mathematical Proceedings of the Cambridge Philosophical Society 26(4), 528–535.