DMRInfer: P-value calculation given Wald statistics.

In ZhenxingGuo0015/TRESS: Toolbox for mRNA epigenetics sequencing analysis

P-value calculation given Wald statistics.

Description

This function calculates p-values for candidate DMRs given Wald statistics, using respectively two-component mixtures of normal, truncated normal and standard normal distributions.

Usage

DMRInfer(stat, nullModel = "standN")

Arguments

stat	A vector containing Wald statistics of all candidate DMRs.
nullModel	A character to specify a method to calculate p-value based on the statistics. It can be "standN", "2mix" and "trunN" for standard normal, two-component mixed gaussian and truncated normal respectively. Defult is "standN".

Details

In addition to standard normal distribution, TRESS provides another two distributions to calculate p-values given Wald statistics, in case that statistics are inflated by potentially underestimated dispersion in data.

One is two-component mixtures of normal distribution, where TRESS assumes that,

T_w \sim p N(0, \sigma_0^2) + (1-p)N(0, \sigma_1^2)

where T_w is Wald statistics, \sigma_0 and \sigma_1 are standard deviations of distribution that T_w follows under null and alternative hypothesis in TRESS_DMRtest. P-values for T_w are calculated using N(0, \hat{\sigma}_0^2).

The other one is truncated normal distribution, where TRESS assumes that,

tT_w\sim p N(0, \sigma^2)

where tT_w = T_w \in [-b, b]. Here, T_w within range [-b, b] is assumed to sampled from null distribuion N(0, \sigma^2). For this truncated normal distribution, TRESS explores different values for boundary b ranging from 1.5 to 2 by step 0.1. TRESS estimates a \hat{\sigma} for each of 6 boundaries. If \hat{\sigma}_{max} - \hat{\sigma}_{min} > 0.5 , TRESS calculates p-values for T_w using N(0, \hat{\sigma}_{min}^2). Otherwise, p-values are obtained using N(0, \hat{\sigma}^2), with \hat{\sigma} estimated under b=2.

Value

This function returns a dataframe containing p-values and Benjamini-Hochberg procedure adjusted p-values.

Examples

### use a randomly generated toy data as an
### illustrate of DMRInfer
set.seed(12345)
p = 0.8
nsites = 10000
flag.TP = rep(NA, nsites)
T_w = rep(NA, nsites)
for (i in seq_len(nsites)) {
  u = runif(1, min = 0, max = 1)
  if(u < p){
    flag.TP[i] = FALSE
    T_w[i] = rnorm(1, 0, sd = 1)
  }else{
    flag.TP[i] = TRUE
    T_w[i] = rnorm(1, 0, sd = 5)
  }
}
res = DMRInfer(stat = T_w, nullModel = "standN")
sum(res$padj < 0.05 & !flag.TP)/sum(res$padj < 0.05)

res = DMRInfer(stat = T_w, nullModel = "2mix")
sum(res$padj < 0.05 & !flag.TP)/sum(res$padj < 0.05)

res = DMRInfer(stat = T_w, nullModel = "trunN")
sum(res$padj < 0.05 & !flag.TP)/sum(res$padj < 0.05)

ZhenxingGuo0015/TRESS documentation built on April 14, 2023, 4:21 p.m.