Browse code

Simplify accessors for batch summaries: Ns, corr, medians, mads

git-svn-id: file:///home/git/hedgehog.fhcrc.org/bioconductor/trunk/madman/Rpacks/crlmm@58646 bc3139a8-67e5-0310-9ffc-ced21a209358

Rob Scharp authored on 01/10/2011 04:48:38
Showing7 changed files

... ...
@@ -1,7 +1,7 @@
1 1
 Package: crlmm
2 2
 Type: Package
3 3
 Title: Genotype Calling (CRLMM) and Copy Number Analysis tool for Affymetrix SNP 5.0 and 6.0 and Illumina arrays.
4
-Version: 1.11.24
4
+Version: 1.11.25
5 5
 Date: 2010-12-10
6 6
 Author: Benilton S Carvalho <Benilton.Carvalho@cancer.org.uk>, Robert Scharpf <rscharpf@jhsph.edu>, Matt Ritchie <mritchie@wehi.edu.au>, Ingo Ruczinski <iruczins@jhsph.edu>, Rafael A Irizarry
7 7
 Maintainer: Benilton S Carvalho <Benilton.Carvalho@cancer.org.uk>, Robert Scharpf <rscharpf@jhsph.edu>, Matt Ritchie <mritchie@wehi.EDU.AU>
... ...
@@ -27,6 +27,7 @@ importMethodsFrom(oligoClasses, allele, calls, "calls<-", confs,
27 27
 		  "confs<-", cnConfidence, "cnConfidence<-", isSnp,
28 28
 		  chromosome, position, A, B,
29 29
 		  "A<-", "B<-", open, close, flags,
30
+		  openff, closeff,
30 31
 		  batchStatistics, "batchStatistics<-", updateObject)
31 32
 importFrom(oligoClasses, chromosome2integer, celfileDate, list.celfiles,
32 33
            copyNumber, initializeBigMatrix, initializeBigVector, isPackageLoaded)
... ...
@@ -46,21 +47,7 @@ importFrom(genefilter, rowSds)
46 47
 importFrom(mvtnorm, dmvnorm)
47 48
 importFrom(ellipse, ellipse)
48 49
 importFrom(ff, ffdf, physical.ff, physical.ffdf, ffrowapply)
49
-## It is important not to import these classes from oligoClasses
50
-## Doing so causes the following errors:
51
-## N.AA(container)[index, ] <- someMatrix
52
-##Error in function (classes, fdef, mtable)  :
53
-##  unable to find an inherited method for function "medianA.AA<-", for signature "CNSet", "ff_matrix"
54
-##importClassesFrom(oligoClasses, ffdf, ff_matrix)
55 50
 importClassesFrom(oligoClasses, ff_matrix, ffdf)
56
-## Important to export these classes
57
-##exportClasses(ff_or_matrix, ff_matrix, ffdf)
58
-##exportClasses(ff_or_matrix)
59
-##---------------------------------------------------------------------------
60
-## lattice imports
61
-##---------------------------------------------------------------------------
62
-##import(panel.number, panel.grid, panel.xyplot, lpoints, lsegments, lrect, ltext)
63
-
64 51
 exportMethods(lines)
65 52
 exportMethods(CA, CB)
66 53
 export(crlmm,
... ...
@@ -1,120 +1,74 @@
1
+
1 2
 setMethod("Ns", signature(object="AssayData"),
2 3
 	  function(object, ...){
3
-		  dotArgs <- names(list(...))
4
-		  missing.i <- !("i" %in% dotArgs)
5
-		  missing.j <- !("j" %in% dotArgs)
6
-		  if(!missing.j) j <- list(...)[["j"]]
7 4
 		  batchnames <- sampleNames(object)
8
-		  if(!missing.j) batchnames <- batchnames[j]
9
-		  if(missing.i & missing.j) stop("Must specify either the rows i or batches j")
10
-		  is.ff <- is(assayDataElement(object, "N.AA"), "ff") | is(assayDataElement(object, "N.AA"), "ffdf")
11
-		  if(is.ff){
12
-			  open(assayDataElement(object, "N.AA"))
13
-			  open(assayDataElement(object, "N.AB"))
14
-			  open(assayDataElement(object, "N.BB"))
15
-		  }
16
-		  N.AA <- as.matrix(assayDataElement(object, "N.AA")[...])
17
-		  N.AB <- as.matrix(assayDataElement(object, "N.AB")[...])
18
-		  N.BB <- as.matrix(assayDataElement(object, "N.BB")[...])
19
-		  if(is.ff){
20
-			  close(assayDataElement(object, "N.AA"))
21
-			  close(assayDataElement(object, "N.AB"))
22
-			  close(assayDataElement(object, "N.BB"))
23
-		  }
5
+		  N.AA <- assayDataElement(object, "N.AA")
6
+		  N.AB <- assayDataElement(object, "N.AB")
7
+		  N.BB <- assayDataElement(object, "N.BB")
24 8
 		  res <- array(NA, dim=c(nrow(N.AA), 3, ncol(N.AA)))
25 9
 		  dimnames(res)[[2]] <- c("AA", "AB", "BB")
26 10
 		  dimnames(res)[[3]] <- batchnames
27
-		  res[, "AA", ] <- N.AA
28
-		  res[, "AB", ] <- N.AB
29
-		  res[, "BB", ] <- N.BB
11
+		  res[, "AA", ] <- N.AA[,]
12
+		  res[, "AB", ] <- N.AB[,]
13
+		  res[, "BB", ] <- N.BB[,]
30 14
 		  return(res)
31 15
 	  })
32 16
 setMethod("corr", signature(object="AssayData"),
33 17
 	  function(object, ...){
34
-		  dotArgs <- names(list(...))
35
-		  missing.i <- !("i" %in% dotArgs)
36
-		  missing.j <- !("j" %in% dotArgs)
37
-		  if(!missing.j) j <- list(...)[["j"]]
38 18
 		  batchnames <- sampleNames(object)
39
-		  if(!missing.j) batchnames <- batchnames[j]
40
-
41
-		  ##if(missing.i & missing.j) stop("Must specify either the rows i or batches j")
42
-		  is.ff <- is(assayDataElement(object, "corrAA"), "ff") | is(assayDataElement(object, "corrAA"), "ffdf")
43
-		  if(is.ff){
44
-			  open(assayDataElement(object, "corrAA"))
45
-			  open(assayDataElement(object, "corrAB"))
46
-			  open(assayDataElement(object, "corrBB"))
47
-		  }
48
-		  corrAA <- as.matrix(assayDataElement(object, "corrAA")[...])
49
-		  corrAB <- as.matrix(assayDataElement(object, "corrAB")[...])
50
-		  corrBB <- as.matrix(assayDataElement(object, "corrBB")[...])
51
-		  if(is.ff){
52
-			  close(assayDataElement(object, "corrAA"))
53
-			  close(assayDataElement(object, "corrAB"))
54
-			  close(assayDataElement(object, "corrBB"))
55
-		  }
19
+		  corrAA <- assayDataElement(object, "corrAA")
20
+		  corrAB <- assayDataElement(object, "corrAB")
21
+		  corrBB <- assayDataElement(object, "corrBB")
56 22
 		  res <- array(NA, dim=c(nrow(corrAA), 3, ncol(corrAA)))
57 23
 		  dimnames(res)[[2]] <- c("AA", "AB", "BB")
58 24
 		  dimnames(res)[[3]] <- batchnames
59
-		  res[, "AA", ] <- corrAA
60
-		  res[, "AB", ] <- corrAB
61
-		  res[, "BB", ] <- corrBB
25
+		  res[, "AA", ] <- corrAA[,]
26
+		  res[, "AB", ] <- corrAB[,]
27
+		  res[, "BB", ] <- corrBB[,]
62 28
 		  return(res)
63 29
 	  })
64 30
 
65 31
 setMethod("medians", signature(object="AssayData"),
66 32
 	  function(object, ...){
67
-		  dotArgs <- names(list(...))
68
-		  missing.i <- !("i" %in% dotArgs)
69
-		  missing.j <- !("j" %in% dotArgs)
70
-		  if(!missing.j) j <- list(...)[["j"]]
71 33
 		  batchnames <- sampleNames(object)
72
-		  if(!missing.j) batchnames <- batchnames[j]
73
-		  if(missing.i & missing.j) stop("Must specify either the rows i or batches j")
74
-		  medianA.AA <- as.matrix(assayDataElement(object, "medianA.AA")[...])
75
-		  medianA.AB <- as.matrix(assayDataElement(object, "medianA.AB")[...])
76
-		  medianA.BB <- as.matrix(assayDataElement(object, "medianA.BB")[...])
77
-		  medianB.AA <- as.matrix(assayDataElement(object, "medianB.AA")[...])
78
-		  medianB.AB <- as.matrix(assayDataElement(object, "medianB.AB")[...])
79
-		  medianB.BB <- as.matrix(assayDataElement(object, "medianB.BB")[...])
34
+		  medianA.AA <- assayDataElement(object, "medianA.AA")
35
+		  medianA.AB <- assayDataElement(object, "medianA.AB")
36
+		  medianA.BB <- assayDataElement(object, "medianA.BB")
37
+		  medianB.AA <- assayDataElement(object, "medianB.AA")
38
+		  medianB.AB <- assayDataElement(object, "medianB.AB")
39
+		  medianB.BB <- assayDataElement(object, "medianB.BB")
80 40
 		  res <- array(NA, dim=c(nrow(medianA.AA), 2, 3, ncol(medianA.AA)))
81 41
 		  dimnames(res)[[2]] <- c("A", "B")
82 42
 		  dimnames(res)[[3]] <- c("AA", "AB", "BB")
83 43
 		  dimnames(res)[[4]] <- batchnames
84
-		  res[, "A", "AA", ] <- medianA.AA
85
-		  res[, "A", "AB", ] <- medianA.AB
86
-		  res[, "A", "BB", ] <- medianA.BB
87
-		  res[, "B", "AA", ] <- medianB.AA
88
-		  res[, "B", "AB", ] <- medianB.AB
89
-		  res[, "B", "BB", ] <- medianB.BB
44
+		  res[, "A", "AA", ] <- medianA.AA[,]
45
+		  res[, "A", "AB", ] <- medianA.AB[,]
46
+		  res[, "A", "BB", ] <- medianA.BB[,]
47
+		  res[, "B", "AA", ] <- medianB.AA[,]
48
+		  res[, "B", "AB", ] <- medianB.AB[,]
49
+		  res[, "B", "BB", ] <- medianB.BB[,]
90 50
 		  return(res)
91 51
 })
92 52
 
93 53
 setMethod("mads", signature(object="AssayData"),
94 54
 	  function(object, ...){
95
-		  dotArgs <- names(list(...))
96
-		  missing.i <- !("i" %in% dotArgs)
97
-		  missing.j <- !("j" %in% dotArgs)
98 55
 		  batchnames <- sampleNames(object)
99
-		  if(!missing.j) j <- list(...)[["j"]]
100
-		  if(!missing.j) batchnames <- batchnames[j]
101
-		  if(missing.i & missing.j) stop("Must specify either the rows i or batches j")
102
-		  madA.AA <- as.matrix(assayDataElement(object, "madA.AA")[...])
103
-		  madA.AB <- as.matrix(assayDataElement(object, "madA.AB")[...])
104
-		  madA.BB <- as.matrix(assayDataElement(object, "madA.BB")[...])
105
-		  madB.AA <- as.matrix(assayDataElement(object, "madB.AA")[...])
106
-		  madB.AB <- as.matrix(assayDataElement(object, "madB.AB")[...])
107
-		  madB.BB <- as.matrix(assayDataElement(object, "madB.BB")[...])
56
+		  madA.AA <- assayDataElement(object, "madA.AA")
57
+		  madA.AB <- assayDataElement(object, "madA.AB")
58
+		  madA.BB <- assayDataElement(object, "madA.BB")
59
+		  madB.AA <- assayDataElement(object, "madB.AA")
60
+		  madB.AB <- assayDataElement(object, "madB.AB")
61
+		  madB.BB <- assayDataElement(object, "madB.BB")
108 62
 		  res <- array(NA, dim=c(nrow(madA.AA), 2, 3, ncol(madA.AA)))
109 63
 		  dimnames(res)[[2]] <- c("A", "B")
110 64
 		  dimnames(res)[[3]] <- c("AA", "AB", "BB")
111 65
 		  dimnames(res)[[4]] <- batchnames
112
-		  res[, "A", "AA", ] <- madA.AA
113
-		  res[, "A", "AB", ] <- madA.AB
114
-		  res[, "A", "BB", ] <- madA.BB
115
-		  res[, "B", "AA", ] <- madB.AA
116
-		  res[, "B", "AB", ] <- madB.AB
117
-		  res[, "B", "BB", ] <- madB.BB
66
+		  res[, "A", "AA", ] <- madA.AA[,]
67
+		  res[, "A", "AB", ] <- madA.AB[,]
68
+		  res[, "A", "BB", ] <- madA.BB[,]
69
+		  res[, "B", "AA", ] <- madB.AA[,]
70
+		  res[, "B", "AB", ] <- madB.AB[,]
71
+		  res[, "B", "BB", ] <- madB.BB[,]
118 72
 		  return(res)
119 73
 })
120 74
 
... ...
@@ -122,25 +76,19 @@ setMethod("mads", signature(object="AssayData"),
122 76
 
123 77
 setMethod("tau2", signature(object="AssayData"),
124 78
 	  function(object, ...){
125
-		  dotArgs <- names(list(...))
126
-		  missing.i <- !("i" %in% dotArgs)
127
-		  missing.j <- !("j" %in% dotArgs)
128 79
 		  batchnames <- sampleNames(object)
129
-		  if(!missing.j) j <- list(...)[["j"]]
130
-		  if(!missing.j) batchnames <- batchnames[j]
131
-		  if(missing.i & missing.j) stop("Must specify either the rows i or batches j")
132
-		  tau2A.AA <- as.matrix(assayDataElement(object, "tau2A.AA")[...])
133
-		  tau2A.BB <- as.matrix(assayDataElement(object, "tau2A.BB")[...])
134
-		  tau2B.AA <- as.matrix(assayDataElement(object, "tau2B.AA")[...])
135
-		  tau2B.BB <- as.matrix(assayDataElement(object, "tau2B.BB")[...])
80
+		  tau2A.AA <- assayDataElement(object, "tau2A.AA")
81
+		  tau2A.BB <- assayDataElement(object, "tau2A.BB")
82
+		  tau2B.AA <- assayDataElement(object, "tau2B.AA")
83
+		  tau2B.BB <- assayDataElement(object, "tau2B.BB")
136 84
 		  res <- array(NA, dim=c(nrow(tau2A.AA), 2, 2, ncol(tau2A.AA)))
137 85
 		  dimnames(res)[[2]] <- c("A", "B")
138 86
 		  dimnames(res)[[3]] <- c("AA", "BB")
139 87
 		  dimnames(res)[[4]] <- batchnames
140
-		  res[, "A", "AA", ] <- tau2A.AA
141
-		  res[, "A", "BB", ] <- tau2A.BB
142
-		  res[, "B", "AA", ] <- tau2B.AA
143
-		  res[, "B", "BB", ] <- tau2B.BB
88
+		  res[, "A", "AA", ] <- tau2A.AA[,]
89
+		  res[, "A", "BB", ] <- tau2A.BB[,]
90
+		  res[, "B", "AA", ] <- tau2B.AA[,]
91
+		  res[, "B", "BB", ] <- tau2B.BB[,]
144 92
 		  return(res)
145 93
 })
146 94
 
... ...
@@ -62,9 +62,9 @@ rawCopynumber(object,...)
62 62
 \value{
63 63
 
64 64
 	nu[A/B] and phi[A/B] return matrices of the intercept and
65
-	slope coefficients, respectively.  
65
+	slope coefficients, respectively.
66 66
 
67
-	CA and CB return matrices of allele-specific copy number.  
67
+	CA and CB return matrices of allele-specific copy number.
68 68
 
69 69
 	totalCopynumber (or rawCopynumber) returns a matrix of CA+CB.
70 70
 
... ...
@@ -79,35 +79,35 @@ rawCopynumber(object,...)
79 79
 \examples{
80 80
 ## Version 1.6* of crlmm used CNSetLM objects.
81 81
 data(sample.CNSet)
82
-all(isCurrent(cnSet)) ## is the cnSet object current?
82
+all(isCurrent(sample.CNSet)) ## is the cnSet object current?
83 83
 
84 84
 ## --------------------------------------------------
85 85
 ## calculating allele-specific copy number
86 86
 ## --------------------------------------------------
87 87
 ## copy number for allele A, first 5 markers, first 2 samples
88
-(ca <- CA(cnSet, i=1:5, j=1:2))
88
+(ca <- CA(sample.CNSet, i=1:5, j=1:2))
89 89
 ## copy number for allele B, first 5 markers, first 2 samples
90
-(cb <- CB(cnSet, i=1:5, j=1:2))
90
+(cb <- CB(sample.CNSet, i=1:5, j=1:2))
91 91
 ## total copy number for first 5 markers, first 2 samples
92 92
 (cn1 <- ca+cb)
93 93
 
94 94
 ## total copy number at first 5 nonpolymorphic loci
95
-index <- which(!isSnp(cnSet))[1:5]
96
-cn2 <- CA(cnSet, i=index, j=1:2)
95
+index <- which(!isSnp(sample.CNSet))[1:5]
96
+cn2 <- CA(sample.CNSet, i=index, j=1:2)
97 97
 ## note, cb is NA at nonpolymorphic loci
98
-(cb <- CB(cnSet, i=index, j=1:2))
98
+(cb <- CB(sample.CNSet, i=index, j=1:2))
99 99
 ## note, ca+cb will give NAs at nonpolymorphic loci
100
-CA(cnSet, i=index, j=1:2) + cb
100
+CA(sample.CNSet, i=index, j=1:2) + cb
101 101
 ## A shortcut for total copy number
102
-cn3 <- totalCopynumber(cnSet, i=1:5, j=1:2)
102
+cn3 <- totalCopynumber(sample.CNSet, i=1:5, j=1:2)
103 103
 all.equal(cn3, cn1)
104
-cn4 <- totalCopynumber(cnSet, i=index, j=1:2)
104
+cn4 <- totalCopynumber(sample.CNSet, i=index, j=1:2)
105 105
 all.equal(cn4, cn2)
106 106
 
107 107
 ## markers 1-5, all samples
108
-cn5 <- totalCopynumber(cnSet, i=1:5)
108
+cn5 <- totalCopynumber(sample.CNSet, i=1:5)
109 109
 ## all markers, samples 1-5
110
-cn6 <- totalCopynumber(cnSet, j=1:5)
110
+cn6 <- totalCopynumber(sample.CNSet, j=1:5)
111 111
 
112 112
 ## NOTE: subsetting the object before extracting copy number
113 113
 ##       can be very inefficient when the data set is very large,
... ...
@@ -116,7 +116,7 @@ cn6 <- totalCopynumber(cnSet, j=1:5)
116 116
 ##       and all of the elements in the LinearModelParameter slot
117 117
 \dontrun{
118 118
         ## do not do the following
119
-	cn <- CA(cnSet[1:5, ], "A")
119
+	cn <- CA(sample.CNSet[1:5, ], "A")
120 120
 }
121 121
 }
122 122
 \keyword{manip}
... ...
@@ -13,8 +13,6 @@ genotype(filenames, cdfName, batch, mixtureSampleSize = 10^5, eps =0.1,
13 13
          verbose = TRUE, seed = 1, sns, probs = rep(1/3, 3),
14 14
          DF = 6, SNRMin = 5, recallMin = 10, recallRegMin = 1000,
15 15
          gender = NULL, returnParams = TRUE, badSNP = 0.7)
16
-
17
-
18 16
 }
19 17
 \arguments{
20 18
   \item{filenames}{ complete path to CEL files}
... ...
@@ -34,7 +32,8 @@ genotype(filenames, cdfName, batch, mixtureSampleSize = 10^5, eps =0.1,
34 32
   \item{recallRegMin}{Minimum number of SNP's for regression.}
35 33
   \item{gender}{  integer vector (  male = 1, female =2 ) or missing,
36 34
   with same length as filenames.  If missing, the gender is predicted.}
37
-  \item{returnParams}{'logical'. Return recalibrated parameters from crlmm.}
35
+  \item{returnParams}{'logical'. Return recalibrated parameters from
36
+    crlmm.}
38 37
   \item{badSNP}{'numeric'. Threshold to flag as bad SNP (affects batchQC)}
39 38
 }
40 39
 
... ...
@@ -91,7 +90,6 @@ if (require(ff) & require(genomewidesnp6Crlmm) & require(hapmapsnp6)){
91 90
   ## the filenames with full path...
92 91
   ## very useful when genotyping samples not in the working directory
93 92
   cels <- list.celfiles(path, full.names=TRUE)
94
-
95 93
   ## Note: one would need at least 10 CEL files for copy number estimation
96 94
   ## To use less RAM, specify a smaller argument to ocProbesets
97 95
   ocProbesets(50e3)
... ...
@@ -108,9 +106,6 @@ if (require(ff) & require(genomewidesnp6Crlmm) & require(hapmapsnp6)){
108 106
   gender <- c("female", "female", "male")
109 107
   gender[gender == "female"] <- 2
110 108
   gender[gender == "male"] <- 1
111
-\dontrun{
112
-  cnSet2 <- (cnSet <- genotype(cels, cdfName="genomewidesnp6", batch=batch, gender=as.integer(gender)))
113
-}
114 109
   dim(cnSet)
115 110
   table(isSnp(cnSet))
116 111
 }
... ...
@@ -1,6 +1,6 @@
1 1
 \name{genotypes}
2 2
 \alias{genotypes}
3
-%- Also NEED an '\alias' for EACH other topic documented here.
3
+
4 4
 \title{
5 5
   The possible genotypes for an integer copy number.
6 6
 }
... ...
@@ -10,15 +10,12 @@
10 10
 \usage{
11 11
 genotypes(copyNumber)
12 12
 }
13
-%- maybe also 'usage' for other objects documented here.
13
+
14 14
 \arguments{
15 15
   \item{copyNumber}{
16
-    Integer (0-4 allowed).
16
+    Integer (0-4 allowed).}
17 17
 }
18
-}
19
-\details{
20 18
 
21
-}
22 19
 \value{
23 20
   Character vector.
24 21
 }
... ...
@@ -27,11 +24,10 @@ genotypes(copyNumber)
27 24
 R. Scharpf
28 25
 }
29 26
 
30
-
31 27
 \examples{
28
+
32 29
 for(i in 0:4) print(genotypes(i))
30
+
33 31
 }
34
-% Add one or more standard keywords, see file 'KEYWORDS' in the
35
-% R documentation directory.
36
-\keyword{manip}
37
-x
32
+
33
+\keyword{manip}
38 34
\ No newline at end of file
... ...
@@ -27,7 +27,7 @@ posteriorProbability(object, predictRegion, copyNumber = 0:4, w)
27 27
 
28 28
   \item{w}{ numeric vector of prior probabilities for each of the copy
29 29
     number states.  Must be the same length as \code{copyNumber} and sum
30
-    to 1.
30
+    to 1.}
31 31
 
32 32
 }
33 33
 \details{