glmTreat {edgeR} | R Documentation |
Conduct genewise statistical tests for a given coefficient or contrast relative to a specified fold-change threshold.
glmTreat(glmfit, coef = ncol(glmfit$design), contrast = NULL, lfc = 0, null = "interval")
glmfit |
a |
coef |
integer or character vector indicating which coefficients of the linear model are to be tested equal to zero. Values must be columns or column names of |
contrast |
numeric vector specifying the contrast of the linear model coefficients to be tested against the log2-fold-change threshold. Length must equal to the number of columns of |
lfc |
numeric scalar specifying the absolute value of the log2-fold change threshold above which differential expression is to be considered. |
null |
character string, choices are |
glmTreat
implements a test for differential expression relative to a minimum required fold-change threshold.
Instead of testing for genes which have log-fold-changes different from zero, it tests whether the log2-fold-change is greater than lfc
in absolute value.
glmTreat
is analogous to the TREAT approach developed by McCarthy and Smyth (2009) for microarrays.
glmTreat
detects whether glmfit
was produced by glmFit
or glmQLFit
.
In the former case, it conducts a modified likelihood ratio test (LRT) against the fold-change threshold.
In the latter case, it conducts a quasi-likelihood (QL) F-test against the threshold.
If lfc=0
, then glmTreat
is equivalent to glmLRT
or glmQLFTest
, depending on whether likelihood or quasi-likelihood is being used.
If there is no shrinkage on log-fold-changes, i.e., fitting glms with prior.count=0
, then unshrunk.logFC
and logFC
are essentially the same. Hence they are merged into one column of logFC
in table
.
Note that glmTreat
constructs test statistics using unshrunk.logFC
rather than logFC
.
glmTreat
with positive lfc
gives larger p-values than would be obtained with lfc=0
.
If null="worst.case"
, then glmTreat
conducts a test closely analogous to the treat
function in the limma package.
This conducts a test if which the null hypothesis puts the true logFC on the boundary of the [-lfc,lfc]
interval closest to the observed logFC.
If null="interval"
, then the null hypotheses assumes an interval of possible values for the true logFC.
This approach is somewhat less conservative.
glmTreat
produces an object of class DGELRT
with the same components as for glmfit
plus the following:
lfc |
absolute value of the specified log2-fold-change threshold. |
table |
data frame with the same rows as |
comparison |
character string describing the coefficient or the contrast being tested. |
The data frame table
contains the following columns:
logFC |
shrunk log2-fold-change of expression between conditions being tested. |
unshrunk.logFC |
unshrunk log2-fold-change of expression between conditions being tested. Exists only when |
logCPM |
average log2-counts per million, the average taken over all libraries. |
PValue |
p-values. |
glmTreat
was previously called treatDGE
.
Yunshun Chen and Gordon Smyth
McCarthy, D. J., and Smyth, G. K. (2009). Testing significance relative to a fold-change threshold is a TREAT. Bioinformatics 25, 765-771. http://bioinformatics.oxfordjournals.org/content/25/6/765
topTags
displays results from glmTreat
.
treat
is the corresponding function in the limma package, designed on normally distributed expression data rather than negative binomial counts.
ngenes <- 100 n1 <- 3 n2 <- 3 nlibs <- n1+n2 mu <- 100 phi <- 0.1 group <- c(rep(1,n1), rep(2,n2)) design <- model.matrix(~as.factor(group)) ### 4-fold change for the first 5 genes i <- 1:5 fc <- 4 mu <- matrix(mu, ngenes, nlibs) mu[i, 1:n1] <- mu[i, 1:n1]*fc counts <- matrix(rnbinom(ngenes*nlibs, mu=mu, size=1/phi), ngenes, nlibs) d <- DGEList(counts=counts,lib.size=rep(1e6, nlibs), group=group) gfit <- glmFit(d, design, dispersion=phi) tr <- glmTreat(gfit, coef=2, lfc=1) topTags(tr)