fdr2d(xdat, grp, test, p0, nperm = 100, nr = 15, seed = NULL, null = NULL,
constrain = TRUE, smooth = 0.2, verb = TRUE, ...)
xdat
xdat
and grp
as the first two
arguments and returns the bivariate test statistics as two-column matrix; by default, two-sample t-statistics and logrithmized standard errors are calculated.test
.tstat
, the test statistic, logse
, the corresponding logarithmized standard error, and fdr.local
, the local false discovery rate. This data frame has the additional class attributes fdr2d.result
and fdr.result
, see Examples. This is the bad old S3 class mechanism employed to provide plot and summary functions.Additional information is provided by a param
attribute, which is a list with the following entries:
p0
was estimated from the data or supplied by the user.fdr1d
. Consequently, many arguments have identical or similar meaning. Specifically for fdr2d
, nr
specifies the number of equidistant breaks defining a two-dimensional grid of cells on which the bivariate test statistics are counted; argument constrain
can be set to ensure that the estimated fdr is decreasing with increasing absolute value of the t-statistic; and argument smooth
specifies the degree of smoothing when estimating the fdr.Note that while fdr2d
might be used for any suitable pair of test statistics, it has only been tested for the default pair, and the smoothing procedure specifically is optimized for this situation.
Note also that the estimation of the proportion p0
directly from the data may be quite unstable and dependant on the degree of smoothing; too heavy smoothing may even lead to estimates greater than 1. It is usually more stable use an estimate of p0
provided by fdr1d
.
Note that fdr1d
can also be used to check the degree of smoothing, see average.fdr
.
plot.fdr2d.result
, summary.fdr.result
, OCshow
, fdr1d
, average.fdr
# We simulate a small example with 5 percent regulated genes and
# a rather large effect size
set.seed(2000)
xdat = matrix(rnorm(50000), nrow=1000)
xdat[1:25, 1:25] = xdat[1:25, 1:25] - 1
xdat[26:50, 1:25] = xdat[26:50, 1:25] + 1
grp = rep(c("Sample A","Sample B"), c(25,25))
# A default run
res2d = fdr2d(xdat, grp)
res2d[1:20,]
# Looking at the results
summary(res2d)
plot(res2d)
res2d[res2d$fdr<0.05, ]
# Extra information
class(res2d)
attr(res2d,"param")
Run the code above in your browser using DataLab