R/calibrate_pvalues.R
calibrate_pvalues.RdRecalibrates a vector of p-values by rescaling each on the 1-df chi-square
scale by a single calibration_factor, mapping a test whose inflation
factor is \(f\) toward \(\lambda = 1\). The same rescaling underlies both
classical genomic control and the LD Score regression intercept correction;
method selects only where \(f\) comes from.
calibrate_pvalues(
p,
method = c("lambda_gc", "ldsc_intercept"),
calibration_factor = NULL,
inflation_only = FALSE
)Numeric vector of p-values.
Calibration method – the source of the factor:
"lambda_gc" (default) or "ldsc_intercept".
Optional numeric scalar, the chi-square divisor
\(f\). For "lambda_gc" it defaults to compute_lambda_gc(p);
for "ldsc_intercept" it is required (the supplied LDSC intercept).
Logical. If TRUE, a (valid) calibration factor
below 1 is replaced by 1, so deflated p-values are left uncorrected –
the field-standard one-sided genomic-control convention (see Details).
Default FALSE: the factor is applied as is, in both directions.
A list with components p (the calibrated vector, same length
and order as the input), calibration_factor (the factor actually
applied; 1 when inflation_only clamped a sub-1 factor), and
method.
On the 1-df chi-square scale the calibration is
$$\chi^2_{\mathrm{cal}} = \chi^2 / f, \qquad
p_{\mathrm{cal}} = \Pr\!\left(\chi^2_1 > \chi^2_{\mathrm{obs}} / f\right),$$
computed in the upper tail for accuracy at small \(p\). The rescaling is
identical across methods; method only selects the source of the factor
\(f\):
"lambda_gc"When calibration_factor is NULL it
defaults to compute_lambda_gc(p) (the genomic-control \(\lambda\)
estimated from p itself). Here calibration_factor is
\(\lambda_{\mathrm{GC}}\), which removes all inflation
(confounding and polygenic).
"ldsc_intercept"calibration_factor is required
(the intercept from ldsc_regression); it is not estimable
from the p-values being calibrated. The intercept removes confounding
while preserving polygenic signal.
A non-finite or non-positive calibration_factor passes p through
unchanged (e.g. "lambda_gc" on fewer than 50 usable p-values, where
compute_lambda_gc returns NA). NA / out-of-\((0,1]\)
entries flow through as NA.
One-sided (inflation-only) convention. Genomic-control correction is
conventionally one-sided: a factor \(f > 1\) (inflation) is corrected, while
\(f < 1\) (deflation) is read as the test's own conservativeness and left
alone – dividing by \(f < 1\) would scale every statistic up and can
manufacture significance. The METAL meta-analysis software has enforced this
since 2008: "a genomic control correction is only applied to studies where the
genomic control parameter lambda > 1.0" (METAL ChangeLog entry of 1 July 2008, where
correcting deflated studies is recorded as a fixed bug; Willer et al. 2010).
Set inflation_only = TRUE to reproduce that convention. The default
FALSE keeps the two-sided rescaling, which applies the factor in both
directions and is useful as a diagnostic of what a naive correction would do.
Bulik-Sullivan B. et al. (2015). LD Score regression distinguishes confounding from polygenicity in GWAS. Nat Genet 47:291.
Willer C.J., Li Y., Abecasis G.R. (2010). METAL: fast and efficient meta-analysis of genomewide association scans. Bioinformatics 26(17):2190-2191. (One-sided convention: its ChangeLog entry of 1 July 2008, https://github.com/statgen/METAL/blob/master/ChangeLog.)
set.seed(1)
chi2 <- rchisq(2000, df = 1) * 1.3 # an inflated null
p <- pchisq(chi2, df = 1, lower.tail = FALSE)
compute_lambda_gc(p) # ~1.3
#> [1] 1.329937
cal <- calibrate_pvalues(p, method = "lambda_gc") # factor estimated from p
compute_lambda_gc(cal$p) # ~1.0
#> [1] 1
# LDSC-intercept correction with an externally estimated factor:
calibrate_pvalues(p, method = "ldsc_intercept", calibration_factor = 1.08)$p[1:3]
#> [1] 0.6257382 0.2478162 0.4804427