setMethod("[", "CrlmmSetList", function(x, i, j, ..., drop = FALSE){
            if (missing(drop)) drop <- FALSE
            if (missing(i) && missing(j))
              {
                 if (length(list(...))!=0)
                   stop("specify genes or samples to subset; use '",
                        substitute(x), "$", names(list(...))[[1]],
                        "' to access phenoData variables")
                 return(x)
               }
            if (!missing(j)){
              f1 <- function(x, j){
                x <- x[, j]
              }
              x <- lapply(x, f1, j)
            }
            if(!missing(i)){
              f2 <- function(x, i){
                x <- x[i, ]
              }
              x <- lapply(x, f2, i)
            }
	as(x, "CrlmmSetList")	
})

setMethod("$", "CrlmmSetList", function(x, name) {
	##if(!(name %in% .parameterNames()[output(x) != 0])){
	if(length(x) != 3){
		stop("'$' operature reserved for accessing parameter names in CopyNumberSet object.  CrlmmSetList must be of length 3")
	}
	j <- grep(name, fvarLabels(x[[3]]))	
	if(length(j) < 1)
		stop(name, " not in fvarLabels of CopyNumberSet object")
	if(length(j) > 1){
		warning("Multiple instances of ", name, " in fvarLabels.  Using the first instance")
		j <- j[1]
	}
	param <- fData(x[[3]])[, j]
	param
})
setMethod(".harmonizeDimnames", "CrlmmSetList", function(object){
	i <- length(object)
	while(i > 1){
		object[[i-1]] <- harmonizeDimnamesTo(object[[i-1]], object[[i]])
		i <- i-1
	}
	object
})

setMethod("A", "CrlmmSetList", function(object) A(object[[1]]))

setMethod("addFeatureAnnotation", "CrlmmSetList", function(object, CHR){
	if(missing(CHR)) stop("Must specificy chromosome")
	cdfName <- annotation(object)
	pkgname <- paste(cdfName, "Crlmm", sep="")	
	path <- system.file("extdata", package=pkgname)
	loader("cnProbes.rda", pkgname=pkgname, envir=.crlmmPkgEnv)
	cnProbes <- get("cnProbes", envir=.crlmmPkgEnv)
	loader("snpProbes.rda", pkgname=pkgname, envir=.crlmmPkgEnv)
	snpProbes <- get("snpProbes", envir=.crlmmPkgEnv)	

	##Feature Data
	snps <- featureNames(object)[snpIndex(object)]
	nps <- featureNames(object)[cnIndex(object)]
	position.snp <- snpProbes[match(snps, rownames(snpProbes)), "position"]
	names(position.snp) <- snps
	position.np <- cnProbes[match(nps, rownames(cnProbes)), "position"]
	names(position.np) <- nps
	
	position <- c(position.snp, position.np)
	position <- position[match(featureNames(object), names(position))]
	stopifnot(identical(names(position), featureNames(object)))
	if(sum(duplicated(names(position))) > 0){
		warning("Removing rows with NA identifiers...")
		##RS: fix this
		I <- which(!is.na(names(position)))
	}  else I <- seq(along=names(position))
	fd <- data.frame(cbind(CHR,
			       position[I]))
	colnames(fd) <- c("chromosome", "position")
	rownames(fd) <- featureNames(object)
	fD <- new("AnnotatedDataFrame",
		  data=fd,
		  varMetadata=data.frame(labelDescription=colnames(fd)))
	return(fD)
})

setMethod("annotation", "CrlmmSetList", function(object) annotation(object[[1]]))
setMethod("B", "CrlmmSetList", function(object) B(object[[1]]))
setMethod("batch", "CrlmmSetList", function(object) batch(object[[3]]))
setMethod("CA", "CrlmmSetList", function(object, ...) CA(object[[3]], ...))
setMethod("CB", "CrlmmSetList", function(object, ...) CB(object[[3]], ...))
setMethod("calls", "CrlmmSetList", function(object) calls(object[[2]]))
setMethod("chromosome", "CrlmmSetList", function(object){
	chr <- NULL
	for(i  in 1:length(object)){
		if(length(fvarLabels(object[[i]])) > 0){
			if("chromosome" %in% fvarLabels(object[[i]])){
				chr <- chromosome(object[[i]])
				break()
			}
		}
	}
	if(is.null(chr)) warning("fvarLabel 'chromosome' not in any element of the CrlmmSetList object")	
	return(chr)
	##chromosome(object[[3]])
	})

setMethod("position", "CrlmmSetList", function(object){
	pos <- NULL
	for(i  in 1:length(object)){
		if(length(fvarLabels(object[[i]])) > 0){
			if("position" %in% fvarLabels(object[[i]])){
				pos <- position(object[[i]])
				break()
			}
		} else next()
	}
	if(is.null(pos)) warning("fvarLabel 'position' not in any element of the CrlmmSetList object")
	return(pos)
	})

setMethod("cnIndex", "CrlmmSetList", function(object, ...) {
	match(cnNames(object[[1]], annotation(object)), featureNames(object))
})
setMethod("combine", signature=signature(x="CrlmmSetList", y="CrlmmSetList"),
          function(x, y, ...){
		  x.abset <- x[[1]]
		  y.abset <- y[[1]]
		  x.snpset <- x[[2]]
		  y.snpset <- y[[2]]
		  abset <- combine(x.abset, y.abset)
		  ##we have hijacked the featureData slot to store parameters.  Biobase will not allow combining our 'feature' data.
		  warning("the featureData is not easily combined...  removing the featureData")
		  ##fd1 <- featureData(x.snpset)
		  ##fd2 <- featureData(y.snpset)
		  featureData(x.snpset) <- annotatedDataFrameFrom(calls(x.snpset), byrow=TRUE)
		  featureData(y.snpset) <- annotatedDataFrameFrom(calls(y.snpset), byrow=TRUE)
		  snpset <- combine(x.snpset, y.snpset)
		  merged <- list(abset, snpset)
		  merged <- as(merged, "CrlmmSetList")
		  merged
	  })
setMethod("confs", "CrlmmSetList", function(object) confs(object[[2]]))
setMethod("copyNumber", "CrlmmSetList", function(object) copyNumber(object[[3]]))
setMethod("dims", "CrlmmSetList", function(object) sapply(object, dim))
setMethod("featureNames", "CrlmmSetList", function(object) featureNames(object[[1]]))
setMethod("ncol", signature(x="CrlmmSetList"), function(x) ncol(x[[1]]))
setMethod("nrow", signature(x="CrlmmSetList"), function(x) nrow(x[[1]]))
setMethod("plot", signature(x="CrlmmSetList"),
	  function(x, y, ...){
		  A <- log2(A(x))
		  B <- log2(B(x))
		  plot(A, B, ...)
	  })
setMethod("points", signature(x="CrlmmSetList"),
	  function(x, y, ...){
		  A <- log2(A(x))
		  B <- log2(B(x))
		  points(A, B, ...)
	  })
setMethod("sampleNames", "CrlmmSetList", function(object) sampleNames(object[[1]]))
setMethod("show", "CrlmmSetList", function(object){
	cat("\n Elements in CrlmmSetList object: \n")
	cat("\n")
	for(i in 1:length(object)){
		cat("class: ", class(object[[i]]), "\n")
		cat("assayData elements: ", ls(assayData(object[[i]])), "\n")
		cat("Dimensions: ", dim(object[[i]]))		
		cat("\n \n")
	}
})
setMethod("snpIndex", "CrlmmSetList", function(object, ...){
	match(snpNames(object[[1]], annotation(object)), featureNames(object))
})
setMethod("splitByChromosome", "CrlmmSetList", function(object, cdfName, outdir){
	path <- system.file("extdata", package=paste(cdfName, "Crlmm", sep=""))	
	load(file.path(path, "snpProbes.rda"))
	load(file.path(path, "cnProbes.rda"))
	k <- grep("chr", colnames(snpProbes))
	if(length(k) < 1) stop("chr or chromosome not in colnames(snpProbes)")
	for(CHR in 1:24){
		cat("Chromosome ", CHR, "\n")
		snps <- rownames(snpProbes)[snpProbes[, k] == CHR]
		cnps <- rownames(cnProbes)[cnProbes[, k] == CHR]
		index <- c(match(snps, featureNames(object)),
			   match(cnps, featureNames(object)))
		index <- index[!is.na(index)]
		crlmmSetList <- object[index, ]
		save(crlmmSetList, file=file.path(outdir, paste("crlmmSetList_", CHR, ".rda", sep="")))
	}
})
setMethod("update", "CrlmmSetList", function(object, ...){
	computeCopynumber(object, ...)
})

setReplaceMethod("CA", signature(object="CrlmmSetList", value="matrix"),
		 function(object, value){
			 CA(object[[3]]) <- value
			 object
		 })
setReplaceMethod("CB", signature(object="CrlmmSetList", value="matrix"),
		 function(object, value){
			 CB(object[[3]]) <- value
			 object
			 })

setReplaceMethod("A", signature(object="CrlmmSetList", value="matrix"),
		 function(object, value){
			 A(object[[1]]) <- value
			 object
		 })
setReplaceMethod("B", signature(object="CrlmmSetList", value="matrix"),
		 function(object, value){
			 B(object[[1]]) <- value
			 object
		 })



setMethod("boxplot", "CrlmmSetList", function(x, ...){
##boxplot.CrlmmSetList <- function(x, ...){
	if(length(x) != 3) stop("elements of list should be of class ABset, SnpSet, and CopyNumberSet, respectively.")
	genotypes <- calls(x)-1
	A1 <- A(x)
	B1 <- B(x)
	Alist <- split(A1, genotypes)
	Alist <- as(rev(Alist), "data.frame")
	Blist <- split(B1, genotypes)
	ylim <- range(unlist(Alist))
	boxplot(Alist, xaxt="n", ylab=expression(I[A]), 
		cex.axis=0.6,
		xlab="",
		ylim=range(unlist(Alist), na.rm=TRUE),	
		border="grey50", xaxs="i", at=0:2, xlim=c(-0.5, 2.5),
		cex.main=0.9, xaxt="n",
		yaxt="n",
		col=cols)
	axis(2, at=pretty(ylim), cex.axis=0.8)
	axis(1, at=0:2, labels=rev(c("2A (AA genotype)", "1A (AB genotype)", "0A (BB genotype)")), cex.axis=0.7)
	##extracts nuA for first batch
	suppressWarnings(nuA <- x$nuA)
	segments(-1, nuA, 
		 0, nuA, lty=2, col="blue")
	suppressWarnings(phiA <- x$phiA)
	##phiA <- fData(x)[["phiA_A"]]
	segments(0, nuA,
		 2.5, nuA+2.5*phiA, lwd=2, col="blue")
	axis(2, at=nuA, labels=expression(hat(nu[A])), cex.axis=0.9)
	text(0, ylim[1], labels=paste("n =", length(Alist[[1]])), cex=0.8)
	text(1, ylim[1], labels=paste("n =", length(Alist[[2]])), cex=0.8)
	text(2, ylim[1], labels=paste("n =", length(Alist[[3]])), cex=0.8)
##	segments(0.5, nuA+0.5*phiA,
##		 1, pretty(ylim, n=5)[2], lty=2, col="blue")
##	segments(1, pretty(ylim, n=5)[2],
##		 1.2, pretty(ylim)[2], lty=2, col="blue")
##	text(1.25, pretty(ylim)[2], pos=4, 
##	     labels=expression(hat(nu[A]) + c[A]*hat(phi[A])))
	legend("topleft", fill=rev(cols), legend=c("AA", "AB", "BB"))##, title="diallelic genotypes")	

	ylim <- range(unlist(Blist))
	boxplot(Blist, xaxt="n", ylab=expression(I[B]), 
		##xlab="diallelic genotypes", 
		cex.axis=0.6, 
		ylim=range(unlist(Blist), na.rm=TRUE),
		border="grey50", xaxs="i", 
		at=0:2, xlim=c(-0.5, 2.5),
		cex.main=0.9, yaxt="n", 
		col=rev(cols))
	axis(2, at=pretty(ylim), cex.axis=0.8)
	axis(1, at=0:2, labels=c("0B (AA genotype)", "1B (AB genotype)", "2B (BB genotype)"), cex.axis=0.7)
	##nuB <- fData(x)[["nuB_A"]]
	suppressWarnings(nuB <- x$nuB)
	segments(-1, nuB, 
		 0, nuB, lty=2, col="blue")
	##phiB <- fData(x[i,])[["phiB_A"]]
	suppressWarnings(phiB <- x$phiB)
	segments(0, nuB,
		 2.5, nuB+2.5*phiB, lwd=2, col="blue")
	axis(2, at=nuB, labels=expression(hat(nu[B])), cex.axis=0.9)
	text(0, ylim[1], labels=paste("n =", length(Blist[[1]])), cex=0.8)
	text(1, ylim[1], labels=paste("n =", length(Blist[[2]])), cex=0.8)
	text(2, ylim[1], labels=paste("n =", length(Blist[[3]])), cex=0.8)
##	segments(0.5, nuB+0.5*phiB,
##		 1, pretty(ylim,n=5)[2], lty=2, col="blue")
##	segments(1, pretty(ylim, n=5)[2],
##		 1.2, pretty(ylim,n=5)[2], lty=2, col="blue")
##	text(1.25, pretty(ylim,n=5)[2], pos=4, labels=expression(hat(nu[B]) + c[B]*hat(phi[B])))	
	mtext(featureNames(x)[i], 3, outer=TRUE, line=1)	
})