| Title: | Controlling bias and inflation in association studies using the empirical null distribution |
|---|---|
| Description: | Bacon can be used to remove inflation and bias often observed in epigenome- and transcriptome-wide association studies. To this end bacon constructs an empirical null distribution using a Gibbs Sampling algorithm by fitting a three-component normal mixture on z-scores. |
| Authors: | Maarten van Iterson [aut, cre], Erik van Zwet [ctb] |
| Maintainer: | Maarten van Iterson <[email protected]> |
| License: | GPL (>= 2) |
| Version: | 1.33.0 |
| Built: | 2024-06-30 02:41:26 UTC |
| Source: | https://github.com/bioc/bacon |
Gibbs Sampler Algorithm to fit a three component normal mixture to z-scores
bacon( teststatistics = NULL, effectsizes = NULL, standarderrors = NULL, niter = 5000L, nburnin = 2000L, nbins = 1000, trim = 0.999, level = 0.05, na.exclude = FALSE, verbose = FALSE, priors = list(sigma = list(alpha = 1.28, beta = 0.36), mu = list(lambda = c(0, 3, -3), tau = c(1000, 100, 100)), epsilon = list(gamma = c(90, 5, 5))), globalSeed = 42, parallelSeed = 42 )bacon( teststatistics = NULL, effectsizes = NULL, standarderrors = NULL, niter = 5000L, nburnin = 2000L, nbins = 1000, trim = 0.999, level = 0.05, na.exclude = FALSE, verbose = FALSE, priors = list(sigma = list(alpha = 1.28, beta = 0.36), mu = list(lambda = c(0, 3, -3), tau = c(1000, 100, 100)), epsilon = list(gamma = c(90, 5, 5))), globalSeed = 42, parallelSeed = 42 )
teststatistics |
numeric vector or matrix of test-statistics |
effectsizes |
numeric vector or matrix of effect-sizes |
standarderrors |
numeric vector or matrix of standard errors |
niter |
number of iterations |
nburnin |
length of the burnin period |
nbins |
default 1000 else bin test-statistics |
trim |
default 0.999 trimming test-statistics |
level |
significance level used to determine prop. null for starting values |
na.exclude |
see ?na.exclude |
verbose |
default FALSE |
priors |
list of parameters for the prior distributions |
globalSeed |
default 42 global seed. If set to NULL, randomization will occur for sequential and parallel bacon calls |
parallelSeed |
default 42 BiocParallel RNGseed. If input statistics are a matrix and globalSeed=NULL, setting parallelSeed=NULL will allow randomization across parallel processes within a bacon call and across separate calls to bacon. |
object of class-Bacon
mvaniterson
Implementation is based on a version from Zhihui Liu https://macsphere.mcmaster.ca/handle/11375/9368
##simulate some test-statistic from a normal mixture ##and run bacon y <- rnormmix(2000, c(0.9, 0, 1, 0, 4, 1)) bc <- bacon(y) ##extract all estimated mixture parameters estimates(bc) ##extract inflation inflation(bc) ##extract bias bias(bc) ##extract bias and inflation corrected test-statistics head(tstat(bc)) ##inspect the Gibbs Sampling output traces(bc) posteriors(bc) fit(bc) ##simulate multiple sets of test-statistic from a normal mixture ##and run bacon y <- matrix(rnormmix(10*2000, c(0.9, 0, 1, 0, 4, 1)), ncol=10) bc <- bacon(y) ##extract all estimated mixture parameters estimates(bc) ##extract only the inflation inflation(bc) ##extract only the bias bias(bc) ##extract bias and inflation corrected P-values head(pval(bc)) ##extract bias and inflation corrected test-statistics head(tstat(bc))##simulate some test-statistic from a normal mixture ##and run bacon y <- rnormmix(2000, c(0.9, 0, 1, 0, 4, 1)) bc <- bacon(y) ##extract all estimated mixture parameters estimates(bc) ##extract inflation inflation(bc) ##extract bias bias(bc) ##extract bias and inflation corrected test-statistics head(tstat(bc)) ##inspect the Gibbs Sampling output traces(bc) posteriors(bc) fit(bc) ##simulate multiple sets of test-statistic from a normal mixture ##and run bacon y <- matrix(rnormmix(10*2000, c(0.9, 0, 1, 0, 4, 1)), ncol=10) bc <- bacon(y) ##extract all estimated mixture parameters estimates(bc) ##extract only the inflation inflation(bc) ##extract only the bias bias(bc) ##extract bias and inflation corrected P-values head(pval(bc)) ##extract bias and inflation corrected test-statistics head(tstat(bc))
An S4 class container for storing Gibbs Sampler input and output
teststatisticsnumeric vector or matrix of test-statistics
effectsizesnumeric vector or matrix of effect-sizes
standarderrorsnumeric vector or matrix of standard errors
tracesarray of Gibbs Sampler traces
estimatesvector or matrix of parameter estimates
priorslist of parameters of for the prior distributions
niternumber of iterations
nburninlength of the burnin period
Method to extract the estimated bias from the 'bacon'-object
bias(object) ## S4 method for signature 'Bacon' bias(object)bias(object) ## S4 method for signature 'Bacon' bias(object)
object |
'bacon'-object |
vector or matrix of inflation
y <- rnormmix(2000, c(0.9, 0, 1, 0, 4, 1)) ##nbins = 100 to speed up the calculations bc <- bacon(y, nbins=100) bias(bc)y <- rnormmix(2000, c(0.9, 0, 1, 0, 4, 1)) ##nbins = 100 to speed up the calculations bc <- bacon(y, nbins=100) bias(bc)
density of a k-component normal mixture
dnormmix(x, theta)dnormmix(x, theta)
x |
x like dnorm(x, ... |
theta |
parameters of the mixture proportion, mean and sd |
details follow
density of a k-component normal mixture
mvaniterson
n <- 2000 theta <- c(0.8, 0, 1, 0, 4, 1) x <- rnormmix(n, theta) hist(x, freq=FALSE, n=100) curve(dnormmix(x, theta), add=TRUE, lwd=2)n <- 2000 theta <- c(0.8, 0, 1, 0, 4, 1) x <- rnormmix(n, theta) hist(x, freq=FALSE, n=100) curve(dnormmix(x, theta), add=TRUE, lwd=2)
Method to extract inflation- and bias-corrected effect-sizes
es(object, corrected = TRUE) ## S4 method for signature 'Bacon' es(object, corrected = TRUE)es(object, corrected = TRUE) ## S4 method for signature 'Bacon' es(object, corrected = TRUE)
object |
'bacon'-object |
corrected |
optional return uncorrected |
vector or matrix of effect-sizes
es <- replicate(6, rnormmix(2000, c(0.9, 0, 1, 0, 4, 1))) se <- replicate(6, 0.8*sqrt(4/rchisq(2000,df=4))) bc <- bacon(NULL, es, se) head(es(bc))es <- replicate(6, rnormmix(2000, c(0.9, 0, 1, 0, 4, 1))) se <- replicate(6, 0.8*sqrt(4/rchisq(2000,df=4))) bc <- bacon(NULL, es, se) head(es(bc))
Method to extract the estimated parameters from the 'bacon'-object
estimates(object) ## S4 method for signature 'Bacon' estimates(object)estimates(object) ## S4 method for signature 'Bacon' estimates(object)
object |
'bacon'-object |
vector or matrix of estimates
y <- rnormmix(2000, c(0.9, 0, 1, 0, 4, 1)) ##nbins = 100 to speed up the calculations bc <- bacon(y, nbins=100) estimates(bc)y <- rnormmix(2000, c(0.9, 0, 1, 0, 4, 1)) ##nbins = 100 to speed up the calculations bc <- bacon(y, nbins=100) estimates(bc)
Method to plot mixture fit
fit(object, index = 1, ...) ## S4 method for signature 'Bacon' fit(object, index, col = "grey75", border = "grey75", ...)fit(object, index = 1, ...) ## S4 method for signature 'Bacon' fit(object, index, col = "grey75", border = "grey75", ...)
object |
'bacon'-object |
index |
if multiple sets of test-statsistics where provided |
... |
additional plotting parameters |
col |
line color default 'grey75' |
border |
border color 'grey75' |
plot of the Gibbs Sampler mixture fit
y <- rnormmix(2000, c(0.9, 0, 1, 0, 4, 1)) ##nbins = 100 to speed up the calculations bc <- bacon(y, nbins=100) fit(bc)y <- rnormmix(2000, c(0.9, 0, 1, 0, 4, 1)) ##nbins = 100 to speed up the calculations bc <- bacon(y, nbins=100) fit(bc)
Method to extract the estimated inflation from the 'bacon'-object
inflation(object) ## S4 method for signature 'Bacon' inflation(object)inflation(object) ## S4 method for signature 'Bacon' inflation(object)
object |
'bacon'-object |
vector or matrix of inflation
y <- rnormmix(2000, c(0.9, 0, 1, 0, 4, 1)) ##nbins = 100 to speed up the calculations bc <- bacon(y, nbins=100) inflation(bc)y <- rnormmix(2000, c(0.9, 0, 1, 0, 4, 1)) ##nbins = 100 to speed up the calculations bc <- bacon(y, nbins=100) inflation(bc)
Perform fixed meta-analysis using inflation and bias corrected effect-sizes and standard errors
meta(object, corrected = TRUE, ...) ## S4 method for signature 'Bacon' meta(object, corrected = TRUE, ...)meta(object, corrected = TRUE, ...) ## S4 method for signature 'Bacon' meta(object, corrected = TRUE, ...)
object |
'bacon'-object |
corrected |
optional return uncorrected |
... |
additional arguments |
TODO maybe add idea's from http://www.netstorm.be/home/meta_analysis#metaAnalysisU
object of class 'bacon' with added fixed-effect meta-analysis test-statistics, effect-sizes and standard-errors
es <- replicate(6, rnormmix(2000, c(0.9, 0, 1, 0, 4, 1))) se <- replicate(6, 0.8*sqrt(4/rchisq(2000,df=4))) bc <- bacon(NULL, es, se) mbc <- meta(bc)es <- replicate(6, rnormmix(2000, c(0.9, 0, 1, 0, 4, 1))) se <- replicate(6, 0.8*sqrt(4/rchisq(2000,df=4))) bc <- bacon(NULL, es, se) mbc <- meta(bc)
simple ggplot2 plotting function for 'bacon'-object
## S4 method for signature 'Bacon' plot(x, y, type = c("hist", "qq"))## S4 method for signature 'Bacon' plot(x, y, type = c("hist", "qq"))
x |
'bacon'-object |
y |
NULL |
type |
hist or qq |
either qq-plot of P-values or histogram of Test-statistics
plot normal mixtures
plotnormmix(x, theta, ...)plotnormmix(x, theta, ...)
x |
vector of test statistics |
theta |
parameters describing the mixture components |
... |
arguments passed to hist |
details follow
return plot with histogram of the data and mixture and individual components
mvaniterson
n <- 2000 theta <- c(0.8, 0, 1, 0, 4, 1) x <- rnormmix(n, theta) plotnormmix(x, theta)n <- 2000 theta <- c(0.8, 0, 1, 0, 4, 1) x <- rnormmix(n, theta) plotnormmix(x, theta)
Method to plot posterior distribution
posteriors( object, thetas = c("sigma.0", "p.0"), index = 1, alphas = c(0.95, 0.9, 0.75), xlab = "", ylab = "", ... ) ## S4 method for signature 'Bacon' posteriors( object, thetas = c("sigma.0", "p.0"), index = 1, alphas = c(0.95, 0.9, 0.75), xlab = "", ylab = "", ... )posteriors( object, thetas = c("sigma.0", "p.0"), index = 1, alphas = c(0.95, 0.9, 0.75), xlab = "", ylab = "", ... ) ## S4 method for signature 'Bacon' posteriors( object, thetas = c("sigma.0", "p.0"), index = 1, alphas = c(0.95, 0.9, 0.75), xlab = "", ylab = "", ... )
object |
'bacon'-object |
thetas |
which thetas to plot |
index |
if multiple sets of test-statsistics where provided |
alphas |
significance level confidence ellipses |
xlab |
optional xlab |
ylab |
optional ylab |
... |
additional plotting parameters |
plot of the Gibbs Sampler posterior probabilities
y <- rnormmix(2000, c(0.9, 0, 1, 0, 4, 1)) ##nbins = 100 to speed up the calculations bc <- bacon(y, nbins=100) posteriors(bc)y <- rnormmix(2000, c(0.9, 0, 1, 0, 4, 1)) ##nbins = 100 to speed up the calculations bc <- bacon(y, nbins=100) posteriors(bc)
Method to extract inflation- and bias-corrected P-values
pval(object, corrected = TRUE) ## S4 method for signature 'Bacon' pval(object, corrected = TRUE)pval(object, corrected = TRUE) ## S4 method for signature 'Bacon' pval(object, corrected = TRUE)
object |
'bacon'-object |
corrected |
optional return uncorrected |
vector or matrix of P-values
y <- rnormmix(2000, c(0.9, 0, 1, 0, 4, 1)) bc <- bacon(y, nbins=100) #nbins = 100 to speed up the calculations head(pval(bc))y <- rnormmix(2000, c(0.9, 0, 1, 0, 4, 1)) bc <- bacon(y, nbins=100) #nbins = 100 to speed up the calculations head(pval(bc))
sample from a normal mixture
rnormmix(n, theta, shuffle = TRUE)rnormmix(n, theta, shuffle = TRUE)
n |
size |
theta |
parameters |
shuffle |
shuffle return vectors or keep nulls and alternative ordered (null, alts) |
details follow
n samples from a normal mixture with parameters theta
mvaniterson
n <- 2000 theta <- c(0.8, 0, 1, 0, 4, 1) x <- rnormmix(n, theta)n <- 2000 theta <- c(0.8, 0, 1, 0, 4, 1) x <- rnormmix(n, theta)
Method to extract inflation- and bias-corrected standard errors
se(object, corrected = TRUE) ## S4 method for signature 'Bacon' se(object, corrected = TRUE)se(object, corrected = TRUE) ## S4 method for signature 'Bacon' se(object, corrected = TRUE)
object |
'bacon'-object |
corrected |
optional return uncorrected |
vector or matrix of standard-errors
es <- replicate(6, rnormmix(2000, c(0.9, 0, 1, 0, 4, 1))) se <- replicate(6, 0.8*sqrt(4/rchisq(2000,df=4))) bc <- bacon(NULL, es, se) head(se(bc))es <- replicate(6, rnormmix(2000, c(0.9, 0, 1, 0, 4, 1))) se <- replicate(6, 0.8*sqrt(4/rchisq(2000,df=4))) bc <- bacon(NULL, es, se) head(se(bc))
Extract top features after meta analysis
topTable( object, number = 10, adjust.method = "bonf", sort.by = c("pval", "eff.size") ) ## S4 method for signature 'Bacon' topTable( object, number = 10, adjust.method = "bonf", sort.by = c("pval", "eff.size") )topTable( object, number = 10, adjust.method = "bonf", sort.by = c("pval", "eff.size") ) ## S4 method for signature 'Bacon' topTable( object, number = 10, adjust.method = "bonf", sort.by = c("pval", "eff.size") )
object |
'bacon'-object |
number |
return specified number of top features, n=-1 return all features |
adjust.method |
P-value multiple testing adjustment method default bonferroni |
sort.by |
order results by pval or eff.size |
table with top features
es <- replicate(6, rnormmix(2000, c(0.9, 0, 1, 0, 4, 1))) se <- replicate(6, 0.8*sqrt(4/rchisq(2000,df=4))) bc <- bacon(NULL, es, se) mbc <- meta(bc) topTable(mbc)es <- replicate(6, rnormmix(2000, c(0.9, 0, 1, 0, 4, 1))) se <- replicate(6, 0.8*sqrt(4/rchisq(2000,df=4))) bc <- bacon(NULL, es, se) mbc <- meta(bc) topTable(mbc)
Method to plot Gibbs sampling traces
traces(object, burnin = TRUE, index = 1) ## S4 method for signature 'Bacon' traces(object, burnin = TRUE, index = 1)traces(object, burnin = TRUE, index = 1) ## S4 method for signature 'Bacon' traces(object, burnin = TRUE, index = 1)
object |
'bacon'-object |
burnin |
include burnin period default true |
index |
if multiple sets of test-statsistics where provided |
plot of the Gibbs Sampler traces
y <- rnormmix(2000, c(0.9, 0, 1, 0, 4, 1)) ##nbins = 100 to speed up the calculations bc <- bacon(y, nbins=100) traces(bc)y <- rnormmix(2000, c(0.9, 0, 1, 0, 4, 1)) ##nbins = 100 to speed up the calculations bc <- bacon(y, nbins=100) traces(bc)
Method to extract inflation- and bias-corrected test-statistics
tstat(object, corrected = TRUE) ## S4 method for signature 'Bacon' tstat(object, corrected = TRUE)tstat(object, corrected = TRUE) ## S4 method for signature 'Bacon' tstat(object, corrected = TRUE)
object |
'bacon'-object |
corrected |
optional return uncorrected |
vector or matrix of test-statistics
y <- rnormmix(2000, c(0.9, 0, 1, 0, 4, 1)) ##nbins = 100 to speed up the calculations bc <- bacon(y, nbins=100) head(tstat(bc))y <- rnormmix(2000, c(0.9, 0, 1, 0, 4, 1)) ##nbins = 100 to speed up the calculations bc <- bacon(y, nbins=100) head(tstat(bc))