topTags {edgeR} | R Documentation |
Extracts the top DE tags in a data frame for a given pair of groups, ranked by p-value or absolute log-fold change.
topTags(object, n=10L, adjust.method="BH", sort.by="PValue", p.value=1)
object |
a |
n |
scalar, number of tags to display/return |
adjust.method |
character string stating the method used to adjust p-values for multiple testing, passed on to |
sort.by |
character string, should the top tags be sorted by p-value ( |
p.value |
cutoff value for adjusted p-values. Only tags with lower p-values are listed. |
an object of class TopTags
containing the following elements for the top n
most differentially expressed tags as determined by sort.by
:
table |
a data frame containing the elements |
adjust.method |
character string stating the method used to adjust p-values for multiple testing. |
comparison |
a vector giving the names of the two groups being compared. |
test |
character string stating the name of the test. |
The dimensions, row names and column names of a TopTags
object are defined by those of table
, see dim.TopTags
or dimnames.TopTags
.
TopTags
objects also have a show
method so that printing produces a compact summary of their contents.
Note that the terms ‘tag’ and ‘gene’ are synonymous here. The function is only named as ‘Tags’ for historical reasons.
Mark Robinson, Davis McCarthy, Gordon Smyth
Robinson MD, Smyth GK (2008). Small-sample estimation of negative binomial dispersion, with applications to SAGE data. Biostatistics 9, 321-332.
Robinson MD, Smyth GK (2007). Moderated statistical tests for assessing differences in tag abundance. Bioinformatics 23, 2881-2887.
Analogous to topTable
in the limma package.
# generate raw counts from NB, create list object y <- matrix(rnbinom(80,size=1,mu=10),nrow=20) d <- DGEList(counts=y,group=rep(1:2,each=2),lib.size=rep(c(1000:1001),2)) rownames(d$counts) <- paste("gene",1:nrow(d$counts),sep=".") # estimate common dispersion and find differences in expression # here we demonstrate the 'exact' methods, but the use of topTags is # the same for a GLM analysis d <- estimateCommonDisp(d) de <- exactTest(d) # look at top 10 topTags(de) # Can specify how many genes to view tp <- topTags(de, n=15) # Here we view top 15 tp # Or order by fold change instead topTags(de,sort.by="logFC")