Browse code

update imputeGender (add args for SNR and SNRMin)

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

Rob Scharp authored on 01/10/2011 04:50:19
Showing4 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.39
4
+Version: 1.11.40
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>
... ...
@@ -85,7 +85,7 @@ crlmmGT <- function(A, B, SNR, mixtureParams, cdfName, row.names=NULL,
85 85
   ##IF gender not provide, we predict
86 86
   if(is.null(gender)){
87 87
 	  if(verbose) message("Determining gender.")
88
-	  gender <- imputeGender(A, B, XIndex, YIndex)
88
+	  gender <- imputeGender(A, B, XIndex, YIndex, SNR, SNRMin)
89 89
   }
90 90
 
91 91
   Indexes <- list(autosomeIndex, XIndex, YIndex)
... ...
@@ -311,7 +311,7 @@ crlmm2 <- function(filenames, row.names=TRUE, col.names=TRUE,
311 311
   return(list2SnpSet(res2, returnParams=returnParams))
312 312
 }
313 313
 
314
-imputeGender <- function(A, B, XIndex, YIndex){
314
+imputeGender <- function(A, B, XIndex, YIndex, SNR, SNRMin){
315 315
 	##XMedian <- apply(log2(A[XIndex,,
316 316
 	##drop=FALSE])+log2(B[XIndex,, drop=FALSE]), 2, median)/2
317 317
 	a <- log2(A[XIndex,,drop=FALSE])
... ...
@@ -50,7 +50,7 @@ crlmmGT2 <- function(A, B, SNR, mixtureParams, cdfName, row.names=NULL,
50 50
 	## FIXME: XIndex may be greater than ocProbesets()
51 51
 	if(is.null(gender)){
52 52
 		if(verbose) message("Determining gender.")
53
-		gender <- imputeGender(A, B, XIndex, YIndex)
53
+		gender <- imputeGender(A, B, XIndex, YIndex, SNR, SNRMin)
54 54
 	}
55 55
 	##
56 56
 	Indexes <- list(autosomeIndex, XIndex, YIndex)
... ...
@@ -89,6 +89,9 @@ if (require(genomewidesnp6Crlmm) & require(hapmapsnp6)){
89 89
   ## very useful when genotyping samples not in the working directory
90 90
   cels <- list.celfiles(path, full.names=TRUE)
91 91
   (crlmmOutput <- crlmm(cels))
92
+  invisible(open(crlmmOutput$gender))
93
+  stopifnot(all(crlmmOutput$gender[] == c(2,2,1)))
94
+  invisible(close(crlmmOutput$gender))
92 95
 
93 96
   ## If gender is known, one should check that the assigned gender is
94 97
   ## correct, or pass the integer coding of gender as an argument to the
... ...
@@ -104,6 +107,12 @@ if (require(genomewidesnp6Crlmm) & require(hapmapsnp6)){
104 107
 	  path <- system.file("celFiles", package="hapmapsnp6")
105 108
 	  cels <- list.celfiles(path, full.names=TRUE)
106 109
 	  gtSet <- genotype(cels, batch=rep(1L,3))
110
+
111
+          gt1 <- calls(crlmmOutput)
112
+          gt2 <- calls(gtSet)[isSnp(gtSet), ]
113
+          gt2 <- gt2[match(rownames(gt2), rownames(gt1)), ]
114
+	  stopifnot(all.equal(gt1, gt2))
115
+
107 116
 	  XIndex <- which(chromosome(gtSet)==23 & isSnp(gtSet))
108 117
 	  YIndex <- which(chromosome(gtSet)==24 & isSnp(gtSet))
109 118
 	  A.X <- log2(A(gtSet)[XIndex,,drop=FALSE])