# Welch t-test p-values are poorly calibrated for $N=2$ samples

I am performing a large number of Welch's t-tests (t-test with unequal variance) on very small sample sizes, often with only two samples per condition. I am finding the p-values are poorly calibrated: the distribution is not uniform on simulated distributions.

Question Is there a known calibration for Welch's t-test for $$N=2$$ samples per condition?

Minimal python code

import numpy as np
import scipy as sp
from matplotlib import pyplot as plt

rng = np.random.default_rng()

X = rng.normal(loc=0, scale=1, size=(2, 1000000))
Y = rng.normal(loc=0, scale=2, size=(2, 1000000))

p = sp.stats.ttest_ind(X, Y, equal_var=False).pvalue
plt.hist(p,
bins=np.linspace(0, 1, 100),
histtype='step',
density=True)
plt.axhline(1)
plt.xlabel('p')
plt.ylabel('PDF of p')

# Statistical test that the distribution is not uniform
ks_stat, ks_p_value = sp.stats.kstest(p, 'uniform')
print("KS Statistic:", ks_stat)
print("KS Test p-value:", ks_p_value)


Here is the output.

KS Statistic: 0.03906382636964301
KS Test p-value: 0.0


• I'm afraid not, considering the difficulty of Behrens-Fisher problem, in such a small sample. But the practical question is why we'd estimate variance using 2 data points? If there is scientific reason to justify a different variance under different conditions, is there additional scientific prior information that can be utilized to analyze the data? Basically, if sample size is small, there is just not much information in it and we have to resort extra knowledge to bring in extra information to assist inference. Commented Jul 17 at 19:34
• Seems to me that with n=2 any procedure that is dealing simultaneously with estimates of two means and two standard deviations is going to perform poorly. I think of it as related to (or the same as) overfitting. Commented Jul 17 at 20:45
• Thanks for all the information! Re scientific reason: I am running an analysis where each sample is one animal, so there are two animals in each condition, and each observation is related to expression of a single gene. The standard procedure is to pool many genes with similar "properties" to build a better estimate of the variance, which is what I am working on (also more mice will be sacrificed in the future). The question arose was to benchmark a Welch t-test against this more complicated procedure, and hence my observation that the p-values have to be calibrated. Commented Jul 19 at 23:58

Question Is there a known calibration for Welch's t-test for N=2 samples per condition?

You can't calibrate this problem because it depends on the unknown ratio of the unequal variances. See Behrens–Fisher problem .

The test is useful anyway because the error/miscalibration is smaller for larger sample sizes.

For n=2, you may consider a random number generator, which gives perfectly uniform p-values. Given that this doesn't hurt the information from the sample much, n=2 is very little, you may wonder about the usefulness of p-values.

• Thanks! I will randomly sample p-values to calibrate the distribution. Commented Jul 22 at 19:14
• To be honest, I don't know the exact details about the current state of the art about Behrens-Fisher problem and I have learnt the passed days that a lot of explanations of Behrens-Fisher distribution are not the same as what Fisher had in mind. For example, the R function bfTest (mentioned in the other answer) is conditioning on the observed ratio of variances, which seems to be not what Fisher had in mind (the alternative, computing the distribution with variable ratio and integrating over it is however a bit more difficult). Commented Jul 22 at 19:21
• Anyway, randomly sampling p-values instead..., I think you are doing fine. Commented Jul 22 at 19:22

(Nice observation!)

The Satterthwaite degrees-of-freedom calculation is only ever claimed to be an approximation, so what you observe is not totally surprising. The behavior seems to vanish as the sample size gets larger. Even with sizes of $$5$$ and $$5$$, the distribution looks uniform to me. If you bump it up to $$8$$ and $$8$$, the KS test p-value is high at $$0.139$$.

A recommendation I have when you look at distributions of p-values to check for uniformity is to look at the empirical CDF.

import numpy as np
import scipy as sp
from matplotlib import pyplot as plt
from statsmodels.distributions.empirical_distribution import ECDF
np.random.seed(2024)

rng = np.random.default_rng()

X = rng.normal(loc=0, scale=1, size=(5, 1000000))
Y = rng.normal(loc=0, scale=2, size=(5, 1000000))

p = sp.stats.ttest_ind(X, Y, equal_var=False).pvalue

ecdf = ECDF(p)(p)
plt.plot([0, 1], [0, 1], c = 'black', label = "U(0, 1)")
plt.scatter(p, ecdf, label = "Empirical CDF")
plt.show()
plt.close()


This makes it easy (at least for me) to see that approximate $$U(0, 1)$$ distribution.

If you really want to do this (which might not necessarily be a good idea), there are solutions to the Behrens-Fisher problem that guarantee the type I error is at most the nominal (typically 5%) value.

The historical references are Fisher 1935, Robinson 1976.

An R implementation is for example the function bfTest in the asht package.

You can go one step further, because the exact PDF of the test statistics is known (but unwieldy), see e.g. Kabe 1966. Having the PDF, you can obtain a one-sided p-value by numerically solving to obtain the quantile function. Annoying, but possible.

A neat list of other possible approaches (to a generalized version of the problem) is at https://fhernanb.github.io/stests/articles/Behrens-Fisher.html

• The term 'exact' for the Behrens-Fisher solution is a bit tricky. The obtained p-values are not exact, what is exact is the computation method and the guarantee to be above a certain minimal level. Commented Jul 18 at 7:01
• @SextusEmpiricus the literature on this is new to me, but my understanding of Kabe's paper is that they provide the exact distribution of the test statistic under the null, so if you used that, the type 1 error rate should be exact. However the commonly implemented approaches don't use the exact distribution (presumably due to difficulties in computation) and so they end up with consevative inexact p-values. Or am I missing something? Commented Jul 18 at 8:42
• The computation gives a fiducial distribution but that doesn't always express an exact confidence distribution. The article from Kane is about computing exactly the fiducial distribution, but not about the exactness of the confidence interval. Because of the lack of space in a comment, I am working this out here: stats.stackexchange.com/a/651332/164061 (it is not finished yet, I want to add a graphic visualisation of the idea why the aproach goes wrong) Commented Jul 18 at 11:59