rowpAUCs-methods {genefilter}R Documentation

Rowwise ROC and pAUC computation

Description

Methods for fast rowwise computation of ROC curves and (partial) area under the curve (pAUC).

Usage

rowpAUCs(x, fac, p=0.1, jitter=FALSE, caseNames=c("1", "2"))

Arguments

x ExpressionSet, exprSet or numeric matrix. The matrix must not contain NA values. Note that the use of exprSet is deprecated.
fac Factor; if x is an ExpressionSet or exprSet, this may also be a character vector of length 1 with the name of a covariate variable in x. fac must have exactly 2 levels.
p Numeric vector of length 1. Limit in (0,1) to integrate pAUC to.
jitter Add small amount of jitter to the data. This can improve the look of the ROC curve when most of the values in x are actually the same since cutpoints are calculated between data values.
caseNames caseNames

{ The names of the two classification cases.}

Details

Rowwise calculation of Receiver Operating Characteristic (ROC) curves and the corresponding partial area under the curve (pAUC) for a given data matrix or ExpressionSet. The function is implemented in C and thus reasonably fast and memory efficient. Cutpoints are calculated before the first, in between and after the last data value. Note: If the total area under the curve is smaller 0.5 it is flipped for plotting and computation of pAUC.

Value

An object of class rowROC with the calculated specificities and sensitivities for each row and the corresponding pAUCs and AUCs values. See rowROC for details.

Methods

Methods exist for rowPAUCs:
rowPAUCs
signature(x="matrix", fac="factor")
rowPAUCs
signature(x="matrix", fac="numeric")
rowPAUCs
signature(x="exprSet")
rowPAUCs
signature(x="exprSet", fac="character")
rowPAUCs
signature(x="ExpressionSet")
rowPAUCs
signature(x="ExpressionSet", fac="character")

Author(s)

Florian Hahne <f.hahne@dkfz.de>

References

Pepe MS, Longton G, Anderson GL, Schummer M.: Selecting differentially expressed genes from microarray experiments. Biometrics. 2003 Mar;59(1):133-42.

See Also

rocdemo.sca, pAUC, rowROC

Examples


   data(sample.ExpressionSet)

   r1 = rowttests(sample.ExpressionSet, "sex")
   r2 = rowpAUCs(sample.ExpressionSet, "sex", p=0.1)

 if(interactive()) {
   plot(area(r2), r1$statistic, pch=16)
   sel <- which(area(r2, TRUE) > 0.7)
   plot(r2[sel])
 }

 ## this compares performance and output of rowpAUCs to function pAUC in
 ## package ROC 
 if(require(ROC)){

  myRule = function(x) {
    pAUC(rocdemo.sca(truth = as.integer(sample.ExpressionSet$sex)-1 ,
data = x, rule = dxrule.sca), t0 = 0.1)
  }
  nGenes = 100
  cat("computation time for ", nGenes, "genes:\n")
  cat("function pAUC: ")
  print(system.time(r3 <- esApply(sample.ExpressionSet[1:nGenes, ], 1, myRule)))
  cat("function rowpAUCs: ")
  print(system.time(r2 <- rowpAUCs(sample.ExpressionSet[1:nGenes, ], "sex")))
  plot(r3, area(r2), xlab="function pAUC", ylab="function rowpAUCs", main="pAUCs")
 }

[Package genefilter version 1.14.1 Index]