Skip to content

Instantly share code, notes, and snippets.

@kaspermunch
Created March 22, 2025 09:21
Show Gist options
  • Select an option

  • Save kaspermunch/19a56792f904421793fac70ba01ea633 to your computer and use it in GitHub Desktop.

Select an option

Save kaspermunch/19a56792f904421793fac70ba01ea633 to your computer and use it in GitHub Desktop.
Fishers method for combining independent p-values
import itertools
import matplotlib.pyplot as plt
from matplotlib.ticker import FuncFormatter
import seaborn as sns
import pandas as pd
import numpy as np
from scipy.stats import chi2
%config InlineBackend.figure_format = 'retina'
def fisher_method(pvalues):
return 1-chi2.cdf(-2 * np.log(pvalues).sum(), 2 * pvalues.size)
def fisher_method_log10(pvalues):
return -np.log10(fisher_method(pvalues))
cmap = 'viridis'
n = 100
x = np.linspace(1e-10, 1, n)
a = np.array(list(itertools.product(x, x)))
b = np.apply_along_axis(fisher_method, 1, a)
c = b.reshape(n, n)
data = pd.DataFrame(c)
data.index = x
data.columns = x
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(12, 5))
sns.heatmap(data=data, cmap=f'{cmap}', cbar_kws={'pad': 0.02, 'label': 'combined p-value'}, ax=ax1)
ax1.contour(np.arange(.5, data.shape[1]), np.arange(.5, data.shape[0]), data, levels=20, colors='white')
ax1.invert_yaxis()
ax1.xaxis.set_major_formatter(FuncFormatter(lambda x, pos: '{:.3g}'.format(x/100)))
ax1.yaxis.set_major_formatter(FuncFormatter(lambda x, pos: '{:.3g}'.format(x/100)))
ax1.set_xlabel('p-value')
ax1.set_ylabel('p-value')
pmin = -8
n = -pmin*1+1
# x = np.logspace(1e-10, 1, n)
x = np.logspace(-8, 0, n)
a = np.array(list(itertools.product(x, x)))
b = np.apply_along_axis(fisher_method_log10, 1, a)
c = b.reshape(n, n)
data = pd.DataFrame(c)
data.index = np.round(-np.log10(x), 1)
data.columns = np.round(-np.log10(x), 1)
sns.heatmap(data=data, cmap=f'{cmap}_r', cbar_kws={'pad': 0.02, 'label': 'combined -log10(p-value)'}, ax=ax2)
cs = ax2.contour(np.arange(.5, data.shape[1]), np.arange(.5, data.shape[0]), data, levels=20, colors='white')
def fmt(x):
return np.round(x,2)
ax2.clabel(cs, cs.levels[::2], inline=True, fmt=fmt, fontsize=10)
ax2.invert_yaxis()
ax2.set_xlabel('-log10(p-value)')
ax2.set_ylabel('-log10(p-value)')
plt.tight_layout()
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment