{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"$\\newcommand{L}[1]{\\| #1 \\|}\\newcommand{VL}[1]{\\L{ \\vec{#1} }}\\newcommand{R}[1]{\\operatorname{Re}\\,(#1)}\\newcommand{I}[1]{\\operatorname{Im}\\, (#1)}$\n",
"\n",
"## Notes on the Bonferroni threshold\n",
"\n",
"The Bonferroni threshold is a family-wise error threshold. That is, it\n",
"treats a set of tests as one *family*, and the threshold is designed to\n",
"control the probability of detecting *any* positive tests in the family\n",
"(set) of tests, if the null hypothesis is true.\n",
"\n",
"### Family-wise error\n",
"\n",
"The Bonferroni correction uses a result from probability theory to\n",
"estimate the probability of finding *any* p value below a threshold\n",
"$\\theta$, given a set (family) of $n$ p values.\n",
"\n",
"When we have found a threshold $\\theta$ that gives a probability\n",
"$\\le \\alpha$ that *any* p value will be $\\lt \\theta$, then\n",
"the threshold $\\theta$ can be said to control the *family-wise\n",
"error rate* at level $\\alpha$.\n",
"\n",
"### Not the Bonferroni correction\n",
"\n",
"The inequality used for the Bonferroni is harder to explain than a\n",
"simpler but related correction, called the Šidák correction.\n",
"\n",
"We will start with that, and then move on to the Bonferroni correction.\n",
"\n",
"The probability that all $n$ tests are *above* p value threshold\n",
"$\\theta$, *assuming tests are independent*:\n",
"\n",
"$$\n",
"(1 - \\theta)^n\n",
"$$\n",
"\n",
"Chance that one or more p values are $\\le \\theta$:\n",
"\n",
"$$\n",
"1 - (1 - \\theta)^n\n",
"$$\n",
"\n",
"We want a uncorrected p value threshold $\\theta$ such that the\n",
"expression above equals some desired family-wise error (FWE) rate\n",
"$\\alpha_{fwe}$. For example we might want a p value threshold\n",
"$\\theta$ such that there is probability ($\\alpha_{fwe}$) of\n",
"0.05 that there is one or more test with $p \\le \\theta$ in a\n",
"family of $n$ tests, on the null hypothesis:\n",
"\n",
"$$\n",
"\\alpha_{fwe} = 1 - (1 - \\theta)^n\n",
"$$\n",
"\n",
"Solve for $\\theta$:\n",
"\n",
"$$\n",
"\\theta = 1 - (1 - \\alpha_{fwe})^{1 / n}\n",
"$$\n",
"\n",
"So, if we have 10 tests, and we want the threshold $\\theta$ to\n",
"control $\\alpha_{fwe}$ at $0.05$:"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"def sidak_thresh(alpha_fwe, n):\n",
" return 1 - (1 - alpha_fwe)**(1./n)\n",
"\n",
"print(sidak_thresh(0.05, 10))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# The Bonferroni correction\n",
"\n",
"$\\newcommand{\\P}{\\mathbb P}$ The Bonferroni correction uses a\n",
"result from probability theory, called Boole’s inequality. The result is\n",
"by George Boole, of *boolean* fame. Boole’s inequality applies to the\n",
"situation where we have a set of events $A_1, A_2, A_3, \\ldots $, each\n",
"with some probability of occurring ${P}(A_1), {P}(A_2), {P}(A_3) \\ldots\n",
"$. The inequality states that the probability of one or more of these\n",
"events occurring is no greater than the sum of the probabilities of the\n",
"individual events:\n",
"\n",
"$$\n",
"\\P\\biggl(\\bigcup_{i} A_i\\biggr) \\le \\sum_i {\\mathbb P}(A_i).\n",
"$$\n",
"\n",
"You can read the $\\cup$ symbol here as “or” or “union”.\n",
"$\\P\\biggl(\\bigcup_{i} A_i\\biggr)$ is the probability of the\n",
"*union* of all events, and therefore the probability of one or more\n",
"event occurring.\n",
"\n",
"Boole’s inequality is true because:\n",
"\n",
"$$\n",
"\\P(A \\cup B) = P(A) + P(B) - P(A \\cap B)\n",
"$$\n",
"\n",
"where you can read $\\cap$ as “and” or “intersection”. Because\n",
"$P(A \\cap B) \\ge 0$:\n",
"\n",
"$$\n",
"\\P(A \\cup B) \\le P(A) + P(B)\n",
"$$\n",
"\n",
"In our case we have $n$ tests (the family of tests). Each test\n",
"that we label as significant is an event. Therefore the sum of the\n",
"probabilities of all possible events is $n\\theta$.\n",
"${\\mathbb P}\\biggl(\\bigcup_{i} A_i\\biggr)$ is our probability of\n",
"family-wise error $\\alpha_{fwe}$. To get a threshold\n",
"$\\theta$ that controls family-wise error at $\\alpha$, we\n",
"need:\n",
"\n",
"$$\n",
"\\frac{\\alpha_{fwe}}{n} \\le \\theta\n",
"$$\n",
"\n",
"For $n=10$ tests and an $\\alpha_{fwe}$ of 0.05:"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"def bonferroni_thresh(alpha_fwe, n):\n",
" return alpha_fwe / n\n",
"\n",
"print(bonferroni_thresh(0.05, 10))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The Bonferroni correction does not assume the tests are independent.\n",
"\n",
"As we have seen, Boole’s inequality relies on:\n",
"\n",
"$$\n",
"\\P(A \\cup B) = P(A) + P(B) - P(A \\cap B) \\implies \\\\\n",
"\\P(A \\cup B) \\le P(A) + P(B)\n",
"$$\n",
"\n",
"This means that the Bonferroni correction will be conservative (the\n",
"threshold will be too low) when the tests are positively dependent\n",
"($P(A \\cap B) \\gg 0$).\n",
"\n",
"The Bonferroni\n",
"$\\theta_{Bonferroni} = \\alpha_{fwe} \\space / \\space n$ is always\n",
"smaller (more conservative) than the Šidák correction\n",
"$\\theta_{Šidák}$ for $n \\ge 1$, but it is close:"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
"np.set_printoptions(precision=4) # print to 4 decimal places\n",
"n_tests = np.arange(1, 11) # n = 1 through 10\n",
"# The exact threshold for independent p values\n",
"print(sidak_thresh(0.05, n_tests))"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# The Bonferroni threshold for the same alpha, n\n",
"print(bonferroni_thresh(0.05, n_tests))"
]
}
],
"metadata": {},
"nbformat": 4,
"nbformat_minor": 2
}