Fit Null Distribution To Censored Data by Maximum Likelihood

censored.fit fits a null distribution to censored data.

fndr.cutoff finds a suitable cutoff point based on the (approximate) false non-discovery rate (FNDR).

censored.fit(x, cutoff, statistic=c("normal", "correlation", "pvalue", "studentt"))
fndr.cutoff(x, statistic=c("normal", "correlation", "pvalue", "studentt"))

Arguments

x

vector of test statistics.

cutoff

truncation point (this may a single value or a vector).

statistic

type of statistic - normal, correlation, or student t.

Details

As null model truncated normal, truncated student t or a truncated correlation density is assumed. The truncation point is specified by the cutoff parameter. All data points whose absolute value are large than the cutoff point are ignored when fitting the truncated null model via maximum likelihood. The total number of data points is only used to estimate the fraction of null values eta0.

Value

censored.fit returns a matrix whose rows contain the estimated parameters and corresponding errors for each cutoff point.

fndr.cutoff returns a tentative cutoff point.

See also

fdrtool.

Examples

# load "fdrtool" library
library("fdrtool")

# simulate normal data
sd.true = 2.232
n = 5000
z = rnorm(n, sd=sd.true)
censored.fit(z, c(2,3,5), statistic="normal")
#>      cutoff N.cens      eta0     eta0.SE       sd      sd.SE
#> [1,]      2   3150 0.9541581 0.010341081 2.094910 0.14425242
#> [2,]      3   4089 0.9788031 0.006533727 2.158073 0.06574381
#> [3,]      5   4883 1.0000000 0.002192344 2.228397 0.03062565


# simulate contaminated mixture of correlation distribution
r = rcor0(700, kappa=10)
u1 = runif(200, min=-1, max=-0.7)
u2 = runif(200, min=0.7, max=1)
rc = c(r, u1, u2)

censored.fit(r, 0.7, statistic="correlation")
#>      cutoff N.cens eta0     eta0.SE    kappa  kappa.SE
#> [1,]    0.7    689    1 0.004787558 9.746592 0.5734004
censored.fit(rc, 0.7, statistic="correlation")
#>      cutoff N.cens      eta0    eta0.SE    kappa  kappa.SE
#> [1,]    0.7    689 0.6379428 0.01485582 9.746592 0.5734004

# pvalue example
data(pvalues)
co = fndr.cutoff(pvalues, statistic="pvalue")
co
#>    cutoff 
#> 0.6767363 
censored.fit(pvalues, cutoff=co, statistic="pvalue")
#>         cutoff N.cens      eta0    eta0.SE
#> [1,] 0.6767363    569 0.4103924 0.01602273