- new function detectionPValues() to compute detection p-values from
negative control probes.
git-svn-id: file:///home/git/hedgehog.fhcrc.org/bioconductor/trunk/madman/Rpacks/limma@118514 bc3139a8-67e5-0310-9ffc-ced21a209358
... | ... |
@@ -1,6 +1,6 @@ |
1 | 1 |
Package: limma |
2 |
-Version: 3.29.7 |
|
3 |
-Date: 2016-06-12 |
|
2 |
+Version: 3.29.8 |
|
3 |
+Date: 2016-06-13 |
|
4 | 4 |
Title: Linear Models for Microarray Data |
5 | 5 |
Description: Data analysis, linear models and differential expression for microarray data. |
6 | 6 |
Author: Gordon Smyth [cre,aut], Yifang Hu [ctb], Matthew Ritchie [ctb], Jeremy Silver [ctb], James Wettenhall [ctb], Davis McCarthy [ctb], Di Wu [ctb], Wei Shi [ctb], Belinda Phipson [ctb], Aaron Lun [ctb], Natalie Thorne [ctb], Alicia Oshlack [ctb], Carolyn de Graaf [ctb], Yunshun Chen [ctb], Mette Langaas [ctb], Egil Ferkingstad [ctb], Marcus Davy [ctb], Francois Pepin [ctb], Dongseok Choi [ctb] |
65 | 67 |
new file mode 100644 |
... | ... |
@@ -0,0 +1,41 @@ |
1 |
+detectionPValues <- function(x,...) UseMethod("detectionPValues") |
|
2 |
+ |
|
3 |
+detectionPValues.EListRaw <- function(x,status=NULL,...) |
|
4 |
+# Detection p-values from negative controls for EListRaw |
|
5 |
+# Gordon Smyth |
|
6 |
+# Created 13 June 2016 |
|
7 |
+{ |
|
8 |
+ if(is.null(status)) status <- x$genes$Status |
|
9 |
+ detectionPValues(x=x$E, status=status, ...) |
|
10 |
+} |
|
11 |
+ |
|
12 |
+detectionPValues.default <- function(x, status, negctrl="negative",...) |
|
13 |
+# Detection p-values from negative controls |
|
14 |
+# Gordon Smyth |
|
15 |
+# Created 13 June 2016 |
|
16 |
+{ |
|
17 |
+ x <- as.matrix(x) |
|
18 |
+ |
|
19 |
+# Identify negative control rows |
|
20 |
+ if(missing(status) || is.null(status)) stop("need to specify status vector") |
|
21 |
+ isneg <- as.integer(status==negctrl) |
|
22 |
+ notneg <- 1L-isneg |
|
23 |
+ nneg <- sum(isneg) |
|
24 |
+ |
|
25 |
+# Count the number of negative controls greater than each value |
|
26 |
+# cs1 counts number of negative controls >= the observed value |
|
27 |
+# cs2 counts number of negative controls > the observed value |
|
28 |
+# By taking the average of cs1 and cs2, we count ties as 1/2 |
|
29 |
+ DetectionPValue <- x |
|
30 |
+ cs1 <- cs2 <- rep_len(0,length.out=nrow(x)) |
|
31 |
+ for (j in 1:ncol(x)) { |
|
32 |
+ o1 <- order(x[,j],isneg,decreasing=TRUE) |
|
33 |
+ o2 <- order(x[,j],notneg,decreasing=TRUE) |
|
34 |
+ cs1[o1] <- cumsum(isneg[o1]) |
|
35 |
+ cs2[o2] <- cumsum(isneg[o2]) |
|
36 |
+ DetectionPValue1[,j] <- cs1+cs2 |
|
37 |
+ } |
|
38 |
+ |
|
39 |
+# Convert to p-values |
|
40 |
+ DetectionPValue/(2*nneg) |
|
41 |
+} |
... | ... |
@@ -26,6 +26,7 @@ There are also a series of utility functions which read the header information f |
26 | 26 |
and produces an \code{ElistRaw} object. |
27 | 27 |
|
28 | 28 |
\code{\link{read.idat}} reads Illumina files in IDAT format, and produces an \code{EListRaw} object. |
29 |
+\code{\link{detectionPValues}} can be used to add detection p-values. |
|
29 | 30 |
|
30 | 31 |
The function \link{as.MAList} can be used to convert a \code{marrayNorm} object to an \code{MAList} object if the data was read and normalized using the marray and marrayNorm packages. |
31 | 32 |
} |
... | ... |
@@ -7,7 +7,7 @@ |
7 | 7 |
Test whether a set of genes is highly ranked relative to other genes in terms of differential expression, accounting for inter-gene correlation. |
8 | 8 |
} |
9 | 9 |
\usage{ |
10 |
-\S3method{camera}{default}(y, index, design, contrast = ncol(design), weights = NULL, |
|
10 |
+\method{camera}{default}(y, index, design, contrast = ncol(design), weights = NULL, |
|
11 | 11 |
use.ranks = FALSE, allow.neg.cor=FALSE, inter.gene.cor=0.01, trend.var = FALSE, |
12 | 12 |
sort = TRUE, \dots) |
13 | 13 |
interGeneCorrelation(y, design) |
14 | 14 |
new file mode 100644 |
... | ... |
@@ -0,0 +1,58 @@ |
1 |
+\name{detectionPValues} |
|
2 |
+\alias{detectionPValues.default} |
|
3 |
+\alias{detectionPValues.EListRaw} |
|
4 |
+\alias{detectionPValues} |
|
5 |
+ |
|
6 |
+\title{Detection P-Values from Negative Controls} |
|
7 |
+ |
|
8 |
+\description{Compute the proportion of negative controls greater than each observed expression value. |
|
9 |
+Particularly useful for Illumina BeadChips.} |
|
10 |
+ |
|
11 |
+\usage{ |
|
12 |
+\method{detectionPValues}{EListRaw}(x, status = NULL, \dots) |
|
13 |
+\method{detectionPValues}{default}(x, status, negctrl = "negative", \dots) |
|
14 |
+} |
|
15 |
+ |
|
16 |
+\arguments{ |
|
17 |
+ \item{x}{object of class \code{EListRaw} or a numeric \code{matrix} containing raw intensities for regular and control probes from a series of microarrays.} |
|
18 |
+ \item{status}{character vector giving probe types. Defaults to \code{x$genes$Status} if \code{x} is an \code{EListRaw} object.} |
|
19 |
+ \item{negctrl}{character string identifier for negative control probes.} |
|
20 |
+ \item{\dots}{other arguments are not currently used.} |
|
21 |
+ } |
|
22 |
+ |
|
23 |
+\details{ |
|
24 |
+The rows of \code{x} for which \code{status == negctrl} are assumed to correspond to negative control probes. |
|
25 |
+ |
|
26 |
+For each column of \code{x}, the detection p-values are defined as \code{(N.eq/2 + N.gt) / N.neg}, where \code{N.gt} is the number of negative controls with expression greater than the observed value, \code{N.eq} is the number of negative controls with expression equal to the observed value, and \code{N.neg} is the total number of negative controls. |
|
27 |
+} |
|
28 |
+ |
|
29 |
+\value{numeric matrix of same dimensions as \code{x} containing detection p-values.} |
|
30 |
+ |
|
31 |
+\references{ |
|
32 |
+Shi, W, de Graaf, C, Kinkel, S, Achtman, A, Baldwin, T, Schofield, L, Scott, H, Hilton, D, Smyth, GK (2010). |
|
33 |
+Estimating the proportion of microarray probes expressed in an RNA sample. |
|
34 |
+\emph{Nucleic Acids Research} 38, 2168-2176. |
|
35 |
+\url{http://nar.oxfordjournals.org/content/38/7/2168} |
|
36 |
+} |
|
37 |
+ |
|
38 |
+\author{Gordon Smyth} |
|
39 |
+ |
|
40 |
+\seealso{ |
|
41 |
+ An overview of LIMMA functions to read expression data is given in \link{03.ReadingData}. |
|
42 |
+ |
|
43 |
+ \code{\link{read.idat}} reads Illumina BeadChip expression data from binary IDAT files. |
|
44 |
+ |
|
45 |
+ \code{\link{neqc}} performs normexp background correction and quantile normalization aided by control probes. |
|
46 |
+} |
|
47 |
+ |
|
48 |
+\examples{ |
|
49 |
+\dontrun{ |
|
50 |
+# Read Illumina binary IDAT files |
|
51 |
+x <- read.idat(idat, bgx) |
|
52 |
+x$genes$DectionPValue <- detectionPValues(x) |
|
53 |
+y <- neqc(x) |
|
54 |
+} |
|
55 |
+} |
|
56 |
+ |
|
57 |
+\keyword{background correction} |
|
58 |
+\keyword{illumina beadchips} |
... | ... |
@@ -8,21 +8,25 @@ propexpr(x, neg.x=NULL, status=x$genes$Status, labels=c("negative","regular")) |
8 | 8 |
\arguments{ |
9 | 9 |
\item{x}{matrix or similar object containing raw intensities for a set of arrays.} |
10 | 10 |
\item{neg.x}{matrix or similar object containing raw intensities for negative control probes for the same arrays. If \code{NULL}, then negative controls must be provided in \code{x}.} |
11 |
- \item{status}{character vector giving probe types.} |
|
12 |
- \item{labels}{character vector giving probe type identifiers.} |
|
11 |
+ \item{status}{character vector specifying control type of each probe. Only used if \code{neg.x} is \code{NULL}.} |
|
12 |
+ \item{labels}{character vector giving the \code{status} values for negative control probes and regular (non-control) probes respectively. If of length 1, then all probes other than the negative controls are assumed to be regular. Only used if \code{neg.x} is \code{NULL}.} |
|
13 | 13 |
} |
14 |
+ |
|
14 | 15 |
\details{ |
15 |
-This function estimates the proportion of expressed in a microarray by utilizing the negative control probes. |
|
16 |
+This function estimates the overall proportion of probes on each microarray that are correspond to expressed genes using the method of Shi et al (2010). |
|
17 |
+The function is especially useful for Illumina BeadChips arrays, although it can in principle be applied to any platform with good quality negative controls. |
|
18 |
+ |
|
19 |
+The negative controls can be supplied either as rows of \code{x} or as a separate matrix. |
|
20 |
+If supplied as rows of \code{x}, then the negative controls are identified by the \code{status} vector. |
|
21 |
+\code{x} might also include other types of control probes, but these will be ignored in the calculation. |
|
22 |
+ |
|
16 | 23 |
Illumina BeadChip arrays contain 750~1600 negative control probes. |
17 |
-The expression profile of these control probes can be saved to a separate file by the Illumina BeadStudio software when using it to output the expression profile for regular probes. |
|
18 |
-The control probe profile could be re-generated if it was not generated when the regular probe profile was created by BeadStudio. |
|
19 |
-Other microarray platforms can also use this function to estimate the proportion of expressed probes in each array, provided that they have a set of negative control probes. |
|
20 |
- |
|
21 |
-\code{labels} can include one or two probe type identifiers. |
|
22 |
-Its first element should be the identifier for negative control probes (\code{negative} by default). |
|
23 |
-If \code{labels} only contains one identifier, then it will be assumed to contain the identifier for negative control probes. |
|
24 |
-By default, \code{regular} is the identifier for regular probes. |
|
24 |
+If \code{read.idat} is used to read Illumina expression IDAT files, then the control probes will be populated as rows of the output \code{EListRaw} object, and the vector \code{x$genes$Status} will be set to identify control probes. |
|
25 |
+ |
|
26 |
+Alternatively, expression values can be exported from Illumina's GenomeStudio software as tab-delimited text files. |
|
27 |
+In this case, the control probes are usually written to a separate file from the regular probes. |
|
25 | 28 |
} |
29 |
+ |
|
26 | 30 |
\value{ |
27 | 31 |
Numeric vector giving the proportions of expressed probes in each array. |
28 | 32 |
} |
... | ... |
@@ -42,7 +46,16 @@ Description to the control probes in Illumina BeadChips can be found in \code{\l |
42 | 46 |
|
43 | 47 |
\examples{ |
44 | 48 |
\dontrun{ |
45 |
-x <- read.ilmn(files="sample probe profile.txt",ctrlfiles="control probe profile.txt") |
|
49 |
+# Read Illumina binary IDAT files |
|
50 |
+x <- read.idat(idat, bgx) |
|
51 |
+propexpr(x) |
|
52 |
+ |
|
53 |
+# Read text files exported from GenomeStudio |
|
54 |
+x <- read.ilmn(files = "sample probe profile.txt", |
|
55 |
+ ctrlfiles = "control probe profile.txt") |
|
46 | 56 |
propexpr(x) |
47 | 57 |
} |
48 | 58 |
} |
59 |
+ |
|
60 |
+\keyword{background correction} |
|
61 |
+\keyword{illumina beadchips} |
... | ... |
@@ -13,14 +13,14 @@ Rotation gene set testing for linear models. |
13 | 13 |
} |
14 | 14 |
|
15 | 15 |
\usage{ |
16 |
-\S3method{roast}{default}(y, index = NULL, design = NULL, contrast = ncol(design), geneid = NULL, |
|
16 |
+\method{roast}{default}(y, index = NULL, design = NULL, contrast = ncol(design), geneid = NULL, |
|
17 | 17 |
set.statistic = "mean", gene.weights = NULL, var.prior = NULL, df.prior = NULL, |
18 | 18 |
nrot = 999, approx.zscore = TRUE, \dots) |
19 |
-\S3method{mroast}{default}(y, index = NULL, design = NULL, contrast = ncol(design), geneid = NULL, |
|
19 |
+\method{mroast}{default}(y, index = NULL, design = NULL, contrast = ncol(design), geneid = NULL, |
|
20 | 20 |
set.statistic = "mean", gene.weights = NULL, var.prior = NULL, df.prior = NULL, |
21 | 21 |
nrot = 999, approx.zscore = TRUE, adjust.method = "BH", |
22 | 22 |
midp = TRUE, sort = "directional", \dots) |
23 |
-\S3method{fry}{default}(y, index = NULL, design = NULL, contrast = ncol(design), geneid = NULL, |
|
23 |
+\method{fry}{default}(y, index = NULL, design = NULL, contrast = ncol(design), geneid = NULL, |
|
24 | 24 |
standardize = "posterior.sd", sort = "directional", \dots) |
25 | 25 |
} |
26 | 26 |
|
... | ... |
@@ -7,7 +7,7 @@ |
7 | 7 |
Gene set enrichment analysis for linear models using rotation tests (ROtation testing using MEan Ranks). |
8 | 8 |
} |
9 | 9 |
\usage{ |
10 |
-\S3method{romer}{default}(y, index, design = NULL, contrast = ncol(design), |
|
10 |
+\method{romer}{default}(y, index, design = NULL, contrast = ncol(design), |
|
11 | 11 |
array.weights = NULL, block = NULL, correlation, |
12 | 12 |
set.statistic = "mean", nrot = 9999, shrink.resid = TRUE, \dots) |
13 | 13 |
} |