Browse code

Incorporating changes recommended by reviewer

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

Benilton Carvalho authored on 28/03/2009 18:14:19
Showing 8 changed files

... ...
@@ -59,6 +59,8 @@ BioC data packages
59 59
 
60 60
 * fixed documentation to acommodate this change
61 61
 
62
-2009-03-28 B Carvalho - committed version 1.0.63
62
+2009-03-28 B Carvalho - committed version 1.0.64
63 63
 
64 64
 * modified warning message
65
+
66
+* modified crlmm and cnrma to look for data at extdata/
... ...
@@ -1,14 +1,14 @@
1 1
 Package: crlmm
2 2
 Type: Package
3 3
 Title: Genotype Calling (CRLMM) and Copy Number Analysis tool for SNP 5.0 and 6.0 arrays.
4
-Version: 1.0.63
4
+Version: 1.0.64
5 5
 Date: 2008-12-28
6 6
 Author: Rafael A Irizarry, Benilton S Carvalho <bcarvalh@jhsph.edu>, Robert Scharpf <rscharpf@jhsph.edu>
7 7
 Maintainer: Benilton S Carvalho <bcarvalh@jhsph.edu>, Robert Scharpf <rscharpf@jhsph.edu>
8 8
 Description: Faster implementation of CRLMM specific to SNP 5.0 and 6.0 arrays, as well as a copy number tool specific to 6.0.
9 9
 License: Artistic-2.0
10 10
 Imports: affyio, preprocessCore, utils, stats, genefilter, splines, Biobase, mvtnorm
11
-Suggests: hapmapsnp5, genomewidesnp5Crlmm, genomewidesnp6Crlmm (>= 1.0.1), methods, GGdata, snpMatrix
11
+Suggests: hapmapsnp5, genomewidesnp5Crlmm (>= 1.0.2), genomewidesnp6Crlmm (>= 1.0.2), methods, GGdata, snpMatrix
12 12
 Collate: AllClasses.R
13 13
          methods-SnpSet.R
14 14
 	 cnrma-functions.R
... ...
@@ -128,7 +128,8 @@ cnrma <- function(filenames, sns, cdfName, seed=1, verbose=FALSE){
128 128
 	require(pkgname, character.only=TRUE) || stop("Package ", pkgname, " not available")
129 129
 	if (missing(sns)) sns <- basename(filenames)
130 130
 	## Loading data in .crlmmPkgEnv and extracting from there
131
-	data("npProbesFid", package=pkgname, envir=.crlmmPkgEnv)
131
+        loader("npProbesFid.rda", .crlmmPkgEnv, pkgname)
132
+##	data("npProbesFid", package=pkgname, envir=.crlmmPkgEnv)
132 133
 	fid <- getVarInEnv("npProbesFid")
133 134
 	gc()
134 135
 	set.seed(seed)
... ...
@@ -141,7 +142,8 @@ cnrma <- function(filenames, sns, cdfName, seed=1, verbose=FALSE){
141 142
 		if (getRversion() > '2.7.0') pb <- txtProgressBar(min=0, max=length(filenames), style=3)
142 143
 	}
143 144
 	##load reference distribution obtained from hapmap
144
-	data(list="1m_reference_cn", package="genomewidesnp6Crlmm", envir=.crlmmPkgEnv)
145
+        loader("1m_reference_cn.rda", .crlmmPkgEnv, pkgname)
146
+##	data(list="1m_reference_cn", package="genomewidesnp6Crlmm", envir=.crlmmPkgEnv)
145 147
 	reference <- getVarInEnv("reference")
146 148
 	for(i in seq(along=filenames)){
147 149
 		y <- as.matrix(read.celfile(filenames[i], intensity.means.only=TRUE)[["INTENSITY"]][["MEAN"]][fid])
... ...
@@ -186,11 +188,14 @@ goodSnps <- function(phi.thr, envir, fewAA=20, fewBB=20){
186 188
 }
187 189
 
188 190
 instantiateObjects <- function(calls, NP, plate, envir, chrom, A, B,
189
-			       gender, SNRmin=5, SNR){
191
+			       gender, SNRmin=5, SNR,
192
+                               pkgname="genomewidesnp6Crlmm"){
190 193
 	envir[["chrom"]] <- chrom
191 194
 	CHR_INDEX <- paste(chrom, "index", sep="")
192
-	data(list=CHR_INDEX, package="genomewidesnp6Crlmm")
193
-        ## Rob, where does 'index[[1]]' come from?
195
+        fname <- paste(CHR_INDEX, ".rda", sep="")
196
+        loader(fname, .crlmmPkgEnv, pkgname)
197
+        index <- get("index", envir=.crlmmPkgEnv)
198
+##	data(list=CHR_INDEX, package="genomewidesnp6Crlmm")
194 199
 	A <- A[index[[1]], SNR > SNRmin]
195 200
 	B <- B[index[[1]], SNR > SNRmin]
196 201
 	calls <- calls[index[[1]], SNR > SNRmin]
... ...
@@ -388,7 +393,7 @@ computeCopynumber <- function(chrom,
388 393
 	}
389 394
 }
390 395
 
391
-nonpolymorphic <- function(plateIndex, NP, envir, CONF.THR=0.99, DF.PRIOR=50){
396
+nonpolymorphic <- function(plateIndex, NP, envir, CONF.THR=0.99, DF.PRIOR=50, pkgname="genomewidesnp6Crlmm"){
392 397
 	p <- plateIndex
393 398
 	CHR <- envir[["chrom"]]
394 399
 	plate <- envir[["plate"]]
... ...
@@ -446,7 +451,9 @@ nonpolymorphic <- function(plateIndex, NP, envir, CONF.THR=0.99, DF.PRIOR=50){
446 451
 	##calculate R2
447 452
 	if(CHR == 23){
448 453
 		cnvs <- envir[["cnvs"]]
449
-		data(cnProbes, package="genomewidesnp6Crlmm")
454
+                loader("cnProbes", pkgname, .crlmmPkgEnv)
455
+                cnProbes <- get("cnProbes", envir=.crlmmPkgEnv)
456
+##		data(cnProbes, package="genomewidesnp6Crlmm")
450 457
 		cnProbes <- cnProbes[match(cnvs, rownames(cnProbes)), ]
451 458
 		par <- cnProbes[, "position"] < 2709520 | (cnProbes[, "position"] > 154584237 & cnProbes[, "position"] < 154913754)
452 459
 		gender <- envir[["gender"]]
... ...
@@ -70,7 +70,9 @@ crlmmGT <- function(A, B, SNR, mixtureParams, cdfName, row.names=NULL,
70 70
   }
71 71
 
72 72
   if(verbose) message("Loading annotations.")
73
-  data(genotypeStuff, mixtureStuff, package=pkgname, envir=.crlmmPkgEnv)
73
+  loader("genotypeStuff.rda", .crlmmPkgEnv, pkgname)
74
+  loader("mixtureStuff.rda", .crlmmPkgEnv, pkgname)
75
+##  data(genotypeStuff, mixtureStuff, package=pkgname, envir=.crlmmPkgEnv)
74 76
 
75 77
   ## this is toget rid of the 'no visible binding' notes
76 78
   ## variable definitions
... ...
@@ -396,9 +398,12 @@ crlmmGTTNoN <- function(A, B, SNR, mixtureParams, cdfName,
396 398
   
397 399
   NC <- ncol(A)
398 400
   NR <- nrow(A)
401
+  pkgname <- getCrlmmAnnotationName(cdfName)
399 402
   
400 403
   if(verbose) message("Loading annotations.")
401
-  data(genotypeStuff, mixtureStuff, package=pkgname, envir=.crlmmPkgEnv)
404
+  loader("genotypeStuff.rda", .crlmmPkgEnv, pkgname)
405
+  loader("mixtureStuff.rda", .crlmmPkgEnv, pkgname)
406
+##  data(genotypeStuff, mixtureStuff, package=pkgname, envir=.crlmmPkgEnv)
402 407
 
403 408
   ## this is toget rid of the 'no visible binding' notes
404 409
   ## variable definitions
... ...
@@ -559,9 +564,12 @@ crlmmGTNormalNoN <- function(A, B, SNR, mixtureParams, cdfName,
559 564
   
560 565
   NC <- ncol(A)
561 566
   NR <- nrow(A)
562
-  
567
+  pkgname <- getCrlmmAnnotationName(cdfName)
568
+
563 569
   if(verbose) message("Loading annotations.")
564
-  data(genotypeStuff, mixtureStuff, package=pkgname, envir=.crlmmPkgEnv)
570
+  loader("genotypeStuff.rda", .crlmmPkgEnv, pkgname)
571
+  loader("mixtureStuff.rda", .crlmmPkgEnv, pkgname)
572
+##  data(genotypeStuff, mixtureStuff, package=pkgname, envir=.crlmmPkgEnv)
565 573
 
566 574
   ## this is toget rid of the 'no visible binding' notes
567 575
   ## variable definitions
... ...
@@ -736,7 +744,9 @@ crlmmGTnm <- function(A, B, SNR, mixtureParams, cdfName, row.names=NULL,
736 744
   }
737 745
 
738 746
   if(verbose) message("Loading annotations.")
739
-  data(genotypeStuff, mixtureStuff, package=pkgname, envir=.crlmmPkgEnv)
747
+  loader("genotypeStuff.rda", .crlmmPkgEnv, pkgname)
748
+  loader("mixtureStuff.rda", .crlmmPkgEnv, pkgname)
749
+##  data(genotypeStuff, mixtureStuff, package=pkgname, envir=.crlmmPkgEnv)
740 750
 
741 751
   ## this is toget rid of the 'no visible binding' notes
742 752
   ## variable definitions
... ...
@@ -14,7 +14,10 @@ snprma <- function(filenames, mixtureSampleSize=10^5, fitMixture=FALSE, eps=0.1,
14 14
   }
15 15
   
16 16
   if(verbose) message("Loading annotations and mixture model parameters.")
17
-  data(preprocStuff, genotypeStuff, mixtureStuff, package=pkgname, envir=.crlmmPkgEnv)
17
+  loader("preprocStuff.rda", .crlmmPkgEnv, pkgname)
18
+  loader("genotypeStuff.rda", .crlmmPkgEnv, pkgname)
19
+  loader("mixtureStuff.rda", .crlmmPkgEnv, pkgname)
20
+##  data(preprocStuff, genotypeStuff, mixtureStuff, package=pkgname, envir=.crlmmPkgEnv)
18 21
   autosomeIndex <- getVarInEnv("autosomeIndex")
19 22
   pnsa <- getVarInEnv("pnsa")
20 23
   pnsb <- getVarInEnv("pnsb")
... ...
@@ -134,3 +134,13 @@ list2SnpSet <- function(x, returnParams=FALSE){
134 134
         data=fd, varMetadata=fdv),
135 135
       annotation=x[["pkgname"]])
136 136
 }
137
+
138
+loader <- function(theFile, envir, pkgname){
139
+  stopifnot(theFile %in% c("genotypeStuff.rda", "mixtureStuff.rda", "preprocStuff.rda"))
140
+  theFile <- file.path(system.file(package=pkgname),
141
+                       "extdata", theFile)
142
+  if (!file.exists(theFile))
143
+    stop("File", theFile, "does not exist in", pkgname)
144
+  load(theFile, envir=envir)
145
+}
146
+
... ...
@@ -1,9 +1,5 @@
1 1
 \name{crlmm}
2 2
 \alias{crlmm}
3
-\alias{calls}
4
-\alias{calls,SnpSet-method}
5
-\alias{confs}
6
-\alias{confs,SnpSet-method}
7 3
 
8 4
 \title{Genotype oligonucleotide arrays with CRLMM}
9 5
 \description{
... ...
@@ -18,12 +14,9 @@ crlmm(filenames, row.names=TRUE, col.names=TRUE,
18 14
       intensityFile, mixtureSampleSize=10^5,
19 15
       eps=0.1, verbose=TRUE, cdfName, sns, recallMin=10,
20 16
       recallRegMin=1000, returnParams=FALSE, badSNP=0.7)
21
-calls(x)
22
-confs(x)
23 17
 }
24 18
 
25 19
 \arguments{
26
-  \item{x}{'SnpSet' object. Usually the result of crlmm.}
27 20
   \item{filenames}{'character' vector with CEL files to be genotyped.}
28 21
   \item{row.names}{'logical'. Use rownames - SNP names?}
29 22
   \item{col.names}{'logical'. Use colnames - Sample names?}
... ...
@@ -73,9 +66,7 @@ if (require(genomewidesnp5Crlmm) & require(hapmapsnp5)){
73 66
   ## the filenames with full path...
74 67
   ## very useful when genotyping samples not in the working directory
75 68
   cels <- list.celfiles(path, full.names=TRUE)
76
-  crlmmOutput <- crlmm(cels)
77
-  (calls(crlmmOutput)[1:10, 1:2])
78
-  (confs(crlmmOutput)[1:10, 1:2])
69
+  (crlmmOutput <- crlmm(cels))
79 70
 }
80 71
 }
81 72
 \keyword{classif}
82 73
new file mode 100644
... ...
@@ -0,0 +1,38 @@
1
+\name{calls}
2
+\alias{calls}
3
+\alias{calls,SnpSet-method}
4
+\alias{confs}
5
+\alias{confs,SnpSet-method}
6
+\title{Accessors for Calls and Confidences on a SnpSet object}
7
+\description{
8
+  \code{calls} returns the genotype calls. CRLMM stores genotype calls
9
+  as integers (1 - AA; 2 - AB; 3 - BB).
10
+
11
+  \code{confs} returns the confidences associated to the genotype
12
+  calls. THe current implementation of CRLMM stores the confidences as
13
+  integers by using the transformation:
14
+
15
+  conf = round(-1000*log2(1-p)),
16
+
17
+  where 'p' is the posterior probability of the call.
18
+}
19
+\usage{
20
+calls(x)
21
+confs(x)
22
+}
23
+\arguments{
24
+  \item{x}{SnpSet object}
25
+}
26
+\value{
27
+  Matrix of genotype calls or confidences.
28
+  }
29
+\examples{
30
+  set.seed(1)
31
+  theCalls <- matrix(sample(1:3, 20, rep=TRUE), nc=2)
32
+  p <- matrix(runif(20), nc=2)
33
+  theConfs <- round(-1000*log2(1-p))
34
+  obj <- new("SnpSet", call=theCalls, callProbability=theConfs)
35
+  calls(obj)
36
+  confs(obj)
37
+}
38
+\keyword{manip}