[Wavetiling-commits] r18 - in pkg: . R

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Wed Feb 15 17:55:39 CET 2012


Author: ppipeler
Date: 2012-02-15 17:55:39 +0100 (Wed, 15 Feb 2012)
New Revision: 18

Modified:
   pkg/NAMESPACE
   pkg/R/allClasses.R
   pkg/R/allGenerics.R
   pkg/R/helperFunctions.R
   pkg/R/initialize-methods.R
   pkg/R/methods-WaveTilingFeatureSet.R
   pkg/R/show-methods.R
Log:


Modified: pkg/NAMESPACE
===================================================================
--- pkg/NAMESPACE	2011-12-21 09:38:39 UTC (rev 17)
+++ pkg/NAMESPACE	2012-02-15 16:55:39 UTC (rev 18)
@@ -44,23 +44,29 @@
 ## Methods from oligo
 ## importMethodsFrom(oligo, ...)
 
+## Plot
+importFrom(graphics, plot)
 
 ## Exporting
 
 ## Classes
-exportClasses(WaveTilingFeatureSet, mapFilterProbe, genomeInfo, Wfm)
+exportClasses(WaveTilingFeatureSet, mapFilterProbe, genomeInfo, WfmFit, WfmInf, WfmFitTime, WfmFitFactor, WfmFitCircadian, WfmFitCustom, WfmInfCompare, WfmInfEffects )
 
 ## exportMethods(show)
+exportMethods(show)
 
+## exportMethods(plot)
+exportMethods(plot)
+
 ## waveTilingFeatureSet methods
 exportMethods( addPheno, getNoGroups, getGroupNames, getReplics, filterOverlap,                      selectProbesFromTilingFeatureSet, bgCorrQn, makeDesign, wfm.analysis)
 
 ## mapFilterProbe methods
 exportMethods(getFilteredIndices, getChromosome, getPosition, getStrand,                            selectProbesFromFilterOverlap)
 
-## Wfm methods
-exportMethods(getProbePosition, getNoProbes, getBetaMAP, getVarBetaMAP, getSmoothPar,               getVarEps, getGenomeInfo, getChromosome, getStrand, getMinPos,                        getMaxPos, getNoLevels, getWfmMethod, getDesign, getPhenoInfo,                        getDataOrigSpace, getDataWaveletSpace, getWaveletFilter, getKj,                       getPrior, getAlpha, getDelta, getTwoSided, getRescale, getSigProbes,                  getRegions, getGenomicRegions, getFDR, getF, getVarF, getEff,                         getVarEff, plotWfm, getNonAnnotatedRegions, getSigGenes)
+## WfmFit methods
+exportMethods(getProbePosition, getNoProbes, getBetaMAP, getVarBetaMAP, getSmoothPar,               getVarEps, getGenomeInfo, getChromosome, getStrand, getMinPos,                        getMaxPos, getNoLevels, getWfmMethod, getDesign, getPhenoInfo,                        getDataOrigSpace, getDataWaveletSpace, getWaveletFilter, getKj,                       getPrior, getAlpha, getDelta, getTwoSided, getRescale, getSigProbes,                  getRegions, getGenomicRegions, getFDR, getF, getVarF, getEff,                         getVarEff, wfm.inference, getSigGenes, getNonAnnotatedRegions)
 
 ## Other
-export(cel2TilingFeatureSet)
+export(cel2TilingFeatureSet, makeContrasts, makeDesign)
 

Modified: pkg/R/allClasses.R
===================================================================
--- pkg/R/allClasses.R	2011-12-21 09:38:39 UTC (rev 17)
+++ pkg/R/allClasses.R	2012-02-15 16:55:39 UTC (rev 18)
@@ -4,9 +4,42 @@
 #genomeInfo
 setClass("genomeInfo",representation(chromosome="vector",strand="character",minPos="numeric",maxPos="numeric"))
 
-#Wfm
-setClass("Wfm",representation(betaMAP="matrix",varbetaMAP="matrix",smoothPar="matrix",varEps="numeric",dataOrigSpace="matrix",dataWaveletSpace="matrix",design="matrix",phenoData="data.frame",method="character",genome.info="genomeInfo",n.levels="numeric",probePosition="vector",wave.filt="character",Kj="numeric",prior="character",alpha="numeric",delta="numeric",two.sided="numeric",rescale="matrix",sigProbes="list",regions="list",GlocRegions="list",FDR="matrix",CI="array",F="matrix",varF="matrix",eff="matrix",varEff="matrix"))
-
 #WaveTilingFeatureSet
 setClass("WaveTilingFeatureSet",contains="TilingFeatureSet")
 
+#WfmFit
+setClass(Class="WfmFit",
+	representation = representation	(betaMAP="matrix",varbetaMAP="matrix",smoothPar="matrix",varEps="numeric",dataOrigSpace="matrix",dataWaveletSpace="matrix",design.matrix="matrix",phenoData="data.frame",genome.info="genomeInfo",n.levels="numeric",probePosition="vector",wave.filt="character",Kj="numeric",prior="character",F="matrix",varF="matrix",P="numeric",Z="matrix",noGroups="numeric",replics="numeric")
+)
+
+### Factor
+setClass(Class="WfmFitFactor", contains="WfmFit")
+
+### Time
+setClass(Class="WfmFitTime", contains="WfmFit")
+
+## Circadian
+setClass(Class="WfmFitCircadian", contains="WfmFit")
+
+## Custom
+setClass(Class="WfmFitCustom", contains="WfmFit")
+
+#WfmInf
+setClass(Class="WfmInf",
+	representation = representation	(alpha="numeric",delta="numeric",two.sided="numeric",sigProbes="list",regions="list",GlocRegions="list",FDR="matrix",CI="array",eff="matrix",varEff="matrix")
+)
+
+## Compare
+setClass(Class="WfmInfCompare", contains="WfmInf")
+
+## Effects
+setClass(Class="WfmInfEffects", contains="WfmInf")
+
+## Means
+setClass(Class="WfmInfMeans", contains="WfmInf")
+
+## Overall Mean
+setClass(Class="WfmInfOverallMean", contains="WfmInf")
+
+## Custom
+setClass(Class="WfmInfCustom", contains="WfmInf")

Modified: pkg/R/allGenerics.R
===================================================================
--- pkg/R/allGenerics.R	2011-12-21 09:38:39 UTC (rev 17)
+++ pkg/R/allGenerics.R	2012-02-15 16:55:39 UTC (rev 18)
@@ -32,19 +32,23 @@
 )
 
 #fit method
-setGeneric("makeDesign",function(object, method=c("twoGroup","compareGroupsTime","compareGroupsFactor","circadian","twoFactors","meansByGroupTime","meansByGroupFactor","effectsTime"), factor.levels=NULL)
+# setGeneric("makeDesign",function(object, method=c("twoGroup","compareGroupsTime","compareGroupsFactor","circadian","twoFactors","meansByGroupTime","meansByGroupFactor","effectsTime"), factor.levels=NULL)
+# {
+# 	standardGeneric("makeDesign")
+# }
+# )
+
+setGeneric("wfm.analysis",function(object, filter.overlap=NULL, design=c("time","circadian","group","factorial","custom"), n.levels, factor.levels=NULL, chromosome, strand, minPos, maxPos, design.matrix=NULL, var.eps=c("margLik","mad"), prior=c("normal","improper"), eqsmooth=TRUE, max.it=20, wave.filt="haar", skiplevels=NULL, trace=FALSE, save.obs=c("plot","regions","all"))
 {
-	standardGeneric("makeDesign")
+	standardGeneric("wfm.analysis")
 }
 )
 
-setGeneric("wfm.analysis",function(object, filter.overlap=NULL, method=c("twoGroup","compareGroupsTime","compareGroupsFactor","circadian","meansByGroupTime","meansByGroupFactor","effectsTime","twoFactors"), n.levels, chromosome, strand, minPos, maxPos, design.matrix=NULL, var.eps=c("margLik","mad"), prior=c("normal","improper"), eqsmooth=TRUE, max.it=20, wave.filt="haar", skiplevels=NULL, trace=FALSE, save.obs=c("plot","regions","all"), alpha=0.05, nsim=1000, delta=NULL, rescale=NULL, two.sided=NULL, minRunPos=90, minRunProbe=1, factor.levels=NULL)
-{
-	standardGeneric("wfm.analysis")
+setGeneric("wfm.inference",function(object,contrast.matrix=NULL, contrasts=c("compare","means","effects","overallMean"), delta=NULL, two.sided=NULL, minRunPos=90, minRunProbe=1, alpha=0.05, nsim=1000, rescale=NULL) {
+	standardGeneric("wfm.inference")
 }
 )
 
-
 #accessors mapFilterProbe
 setGeneric("getFilteredIndices",function(object)
 {
@@ -161,6 +165,11 @@
 	standardGeneric("getDesign")
 }
 )
+setGeneric("getWfmDesign",function(object)
+{
+	standardGeneric("getWfmDesign")
+}
+)
 setGeneric("getPhenoInfo",function(object)
 {
 	standardGeneric("getPhenoInfo")
@@ -256,19 +265,45 @@
 	standardGeneric("getVarEff")
 }
 )
+setGeneric("getDesignMatrix",function(object)
+{
+	standardGeneric("getDesignMatrix")
+}
+)
 
+setGeneric("getP",function(object)
+{
+	standardGeneric("getP")
+}
+)
+
+setGeneric("getZ",function(object)
+{
+	standardGeneric("getZ")
+}
+)
 setGeneric("plotWfm",function(object, annoFile, minPos, maxPos, trackFeature="exon", overlayFeature=c("gene","transposable_element_gene"), two.strand=TRUE, plotData=TRUE, plotMean=TRUE, tracks=0)
 {
 	standardGeneric("plotWfm")
 }
 )
-setGeneric("getNonAnnotatedRegions",function(object, annoFile)
+
+setGeneric("getSigGenes",function(fit, inf, annoFile)
 {
-	standardGeneric("getNonAnnotatedRegions")
+	standardGeneric("getSigGenes")
 }
 )
-setGeneric("getSigGenes",function(object, annoFile)
+
+setGeneric("getNonAnnotatedRegions",function(fit, inf, annoFile)
 {
-	standardGeneric("getSigGenes")
+	standardGeneric("getNonAnnotatedRegions")
 }
 )
+
+if (!isGeneric("plot")) {
+  setGeneric("plot",function(fit,inf,...){
+	standardGeneric("plot")
+  }
+) 
+}
+

Modified: pkg/R/helperFunctions.R
===================================================================
--- pkg/R/helperFunctions.R	2011-12-21 09:38:39 UTC (rev 17)
+++ pkg/R/helperFunctions.R	2012-02-15 16:55:39 UTC (rev 18)
@@ -284,3 +284,370 @@
 	return(makeRectangleOverlay(start=locations[detectedRegionsSelect[,1]],end =locations[detectedRegionsSelect[,2]],region=region,dp=dp)) 
 }
 
+#### makeDesign
+makeDesign<-function(design=c("time","circadian","group","factorial"),replics, noGroups, factor.levels=NULL) {
+# 	if (method=="twoGroup")
+# 	{
+# 		Xorig <- matrix(0,nrow=dim(pData(object))[1],ncol=2)
+# 		Xorig[,1] <- 1
+# 		Xorig[,2] <- rep(c(1,-1),replics) ## contr.helmert
+# 	}
+	if (design=="time") #time
+	{
+#		Xorig <- matrix(0,nrow=dim(pData(object))[1],ncol=noGroups)
+		Xorig <- matrix(0,nrow=noGroups*replics,ncol=noGroups)
+		orderedFactor <- factor(1:noGroups,ordered=TRUE)
+		desPoly <- lm(rnorm(noGroups)~orderedFactor,x=TRUE)$x
+		Xorig[,1] <- 1
+		Xorig[,2:noGroups] <- apply(desPoly[,2:noGroups],2,rep,replics)
+		
+	}
+	else if (design=="circadian") # circadian
+	{
+#		Xorig <- matrix(0,nrow=dim(pData(object))[1],ncol=3)
+		Xorig <- matrix(0,nrow=noGroups*replics,ncol=3)
+		Xorig[,1] <- 1
+		Xorig[,2] <- rep(sin(seq(0,2*pi-(pi/noGroups),2*pi/noGroups)),replics)
+		Xorig[,3] <- rep(cos(seq(0,2*pi-(pi/noGroups),2*pi/noGroups)),replics)
+	}
+	else  if (design=="group") # group
+	{
+#		Xorig <- matrix(0,nrow=dim(pData(object))[1],ncol=noGroups)
+		Xorig <- matrix(0,nrow=noGroups*replics,ncol=noGroups)
+		desHelmert <- contr.helmert(noGroups)
+		Xorig[,1] <- 1
+		Xorig[,2:noGroups] <- apply(desHelmert[,1:(noGroups-1)],2,rep,replics)
+	}
+	else if (method=="factorial")
+	{
+		if (is.null(factor.levels))
+		{
+			stop("No proper factor levels given.")
+		}
+		Xorig <- matrix(0,nrow=sum(replics),ncol=prod(factor.levels))
+		Xorig[,1] <- 1
+		desTrt1 <- contr.treatment(factor.levels[1])
+		desTrt2 <- contr.treatment(factor.levels[2])
+		#desHelmert1 <- contr.helmert(factor.levels[1])
+		#desHelmert2 <- contr.helmert(factor.levels[2])
+		replicsMat <- matrix(replics,nrow=factor.levels[1],ncol=factor.levels[2],byrow=TRUE)
+		desFac1 <- apply(desTrt1,2,function(x) unlist(sapply(1:nrow(replicsMat),function(y) rep(x[y],sum(replicsMat[y,])))))
+		Xorig[,2:(factor.levels[1])] <- desFac1
+		desFac2 <- apply(desTrt2,2,function(x) unlist(sapply(1:nrow(replicsMat),function(y) rep(x,replicsMat[y,]))))
+		Xorig[,(factor.levels[1]+1):(factor.levels[1]+factor.levels[2]-1)] <- desFac2
+		# include interactions
+		colno <- sum(factor.levels)
+		for (i in 1:ncol(desFac1))
+		{
+			for (j in 1:ncol(desFac2))
+			{
+				Xorig[,colno] <- desFac1[,i]*desFac2[,j]
+				colno <- colno + 1
+			}
+		}
+	}
+	else {
+		stop("No proper design specified!")
+	}	
+	out <- list()
+	X.qr <- qr(Xorig)
+	out$Xorthnorm <- qr.Q(X.qr)
+	out$Xorig <- Xorig
+	return (out);
+}
+
+#### makeContrasts
+makeContrasts <- function(contrasts, nlevels) {
+	if (contrasts=="compare") {  
+		q <- nlevels*(nlevels-1)/2
+		contr <- matrix(0,nrow=q,ncol=nlevels)
+		hlp1 <- rep(2:nlevels,1:(nlevels-1))
+		hlp2 <- unlist(sapply(1:(nlevels-1),function(x) seq(1:x)))
+		for (i in 1:nrow(contr))
+		{
+			contr[i,hlp1[i]] <- 1
+			contr[i,hlp2[i]] <- -1
+		}
+	}
+
+	return (contr);
+}
+
+
+##### plot
+plotWfm<-function(fit,inf,annoFile,minPos,maxPos,trackFeature="exon",overlayFeature=c("gene","transposable_element_gene"),two.strand=TRUE,plotData=TRUE,plotMean=TRUE,tracks=0)
+{
+	if (missing(annoFile)) {stop("Annotation File is missing!!")}
+	Gloc <- getProbePosition(fit)
+	if (missing(minPos)) {minPos<-min(Gloc)}
+	if (missing(maxPos)) {maxPos<-max(Gloc)}
+	chromosome <- getChromosome(fit)
+	strand <- getStrand(fit)
+	selID <- (1:length(Gloc))[Gloc>minPos & Gloc<maxPos]
+	sta <- min(selID)
+	end <- max(selID)
+	minBase <- Gloc[sta]
+	maxBase <- Gloc[end]
+	trackCount <- 1
+	trackInfo <- list()
+	overlayInfo <- list()
+	if (plotData==TRUE)
+	{
+		Y <- getDataOrigSpace(fit)
+		replcs <- as.numeric(table(getPhenoInfo(fit)$group))
+		trackInfo[[trackCount]] <- makeGenericArray(intensity=as.matrix(t(Y[,sta:end])),probeStart=Gloc[sta:end],dp=DisplayPars(color=rep(1:length(replcs),replcs),ylim=range(Y[,sta:end]),pointSize=.3,pch=sequence(replcs)))
+		names(trackInfo)[trackCount] <- "Data"
+		trackCount <- trackCount + 1
+	}
+	if (two.strand==TRUE)
+	{
+		trackInfo[[trackCount]] <- makeNewAnnotationTrack(annoFile=annoFile,chromosome=chromosome,minBase=minBase,maxBase=maxBase,strand="forward",feature=trackFeature,dp=NULL)
+		names(trackInfo)[trackCount] <- "F"
+		overlayInfo[[trackCount]] <- makeNewAnnotationTextOverlay(annoFile=annoFile,chromosome=chromosome,minBase=minBase,maxBase=maxBase,strand="forward",region=c(trackCount,trackCount),feature=overlayFeature,y=0.5)
+		trackCount <- trackCount + 1
+		trackInfo[[trackCount]] <- makeGenomeAxis(add53 = TRUE,add35 = TRUE)
+		trackCount <- trackCount + 1
+		trackInfo[[trackCount]] <- makeNewAnnotationTrack(annoFile=annoFile,chromosome=chromosome,minBase=minBase,maxBase=maxBase,strand="reverse",feature=trackFeature)
+		names(trackInfo)[trackCount] <- "R"
+		overlayInfo[[trackCount]] <- makeNewAnnotationTextOverlay(annoFile=annoFile,chromosome=chromosome,minBase=minBase,maxBase=maxBase,strand="reverse",region=c(trackCount,trackCount),feature=overlayFeature,y=0.5)
+		trackCount <- trackCount + 1
+	} else
+	{
+		trackInfo[[trackCount]] <- makeNewAnnotationTrack(annoFile=annoFile,chromosome=chromosome,minBase=minBase,maxBase=maxBase,strand=strand,feature=trackFeature)
+		if (strand=="forward")
+		{
+			names(trackInfo)[trackCount] <- "F"
+		}
+		if (strand=="reverse")
+		{
+			names(trackInfo)[trackCount] <- "R"
+		}
+		overlayInfo[[trackCount]] <- makeNewAnnotationTextOverlay(annoFile=annoFile,chromosome=chromosome,minBase=minBase,maxBase=maxBase,strand=strand,region=c(trackcount,trackCount),feature=overlayFeature,y=0.5)
+		trackCount <- trackCount + 1
+		gAxis <- makeGenomeAxis(add53 = TRUE,add35 = TRUE)
+		trackCount <- trackCount + 1
+	}
+	effects <- getEff(inf)
+	regions <- getRegions(inf)
+	noGroups <- length(table(getPhenoInfo(fit)$group))
+	if (inherits(inf,"WfmInfMeans")) {
+		plotMean <- FALSE
+	}
+	if (plotMean==TRUE)
+	{
+		trackInfo[[trackCount]] <- makeGenericArray(intensity=as.matrix(effects[1,sta:end]),probeStart=Gloc[sta:end],dp=DisplayPars(color="black",ylim=range(effects[1,sta:end])+c(-0.4,0.4),pointSize=.3,pch=1,lwd=1,type="line"))
+		names(trackInfo)[trackCount] <- "Mean"
+		overlayInfo[[trackCount]] <- makeNewTranscriptRectangleOverlay(sigRegions=as.matrix(data.frame(start(regions[[1]]),end(regions[[1]]))),location=Gloc,start=sta,end=end,region=c(trackCount,trackCount),dp=DisplayPars(color="darkgrey",alpha=.1))
+		trackCount <- trackCount + 1
+	}
+#	if (getWfmFitMethod(object)=="twoGroup" | getWfmFitMethod(object)=="circadian")
+	if ((inherits(fit,"WfmFitFactor") & noGroups==2) | inherits(fit,"WfmFitCircadian"))  
+	{
+		if (tracks==0)
+		{
+			trackInfo[[trackCount]] <- makeGenericArray(intensity=as.matrix(effects[2,sta:end]),probeStart=Gloc[sta:end],dp=DisplayPars(color="black",ylim=c(min(-1.3,range(effects[2,sta:end])[1]),max(1.3,range(effects[2,sta:end])[2]))+c(-0.2,0.2),pointSize=.3,pch=1,lwd=1,type="line"))
+			#if (getWfmFitMethod(object)=="twoGroup")
+			if (inherits(inf,"WfmInfCompare")) 
+			{
+				trackInfo[[trackCount]] <- makeGenericArray(intensity=as.matrix(effects[2,sta:end]),probeStart=Gloc[sta:end],dp=DisplayPars(color="black",ylim=c(min(-1.3,range(effects[2,sta:end])[1]),max(1.3,range(effects[2,sta:end])[2])+c(-0.2,0.2)),pointSize=.3,pch=1,lwd=1,type="line"))
+				names(trackInfo)[trackCount] <- "FC"
+			} else
+			{
+				trackInfo[[trackCount]] <- makeGenericArray(intensity=as.matrix(sqrt(effects[2,sta:end]^2+effects[3,sta:end]^2)),probeStart=Gloc[sta:end],dp=DisplayPars(color="black",ylim=c(min(-1.3,range(sqrt(effects[2,sta:end]^2+effects[3,sta:end]^2))[1]),max(1.2,range(sqrt(effects[2,sta:end]^2+effects[3,sta:end]^2))))+c(-0.2,0.2),pointSize=.3,pch=1,lwd=1,type="line"))
+				names(trackInfo)[trackCount] <- "Ampl"
+			}
+			overlayInfo[[trackCount]] <- makeNewTranscriptRectangleOverlay(sigRegions=as.matrix(data.frame(start(regions[[2]]),end(regions[[2]]))),location=Gloc,start=sta,end=end,region=c(trackCount,trackCount),dp=DisplayPars(color="darkgrey",alpha=.1))
+			hlp <- which(!unlist(lapply(overlayInfo,is.null)))
+			if (two.strand==TRUE)
+			{
+				hlpAnno <- hlp[1:2]
+			} else
+			{
+				hlpAnno <- hlp[1]
+			}
+			hlpAnnoIndex <- unlist(lapply(hlpAnno,function(x) length(overlayInfo[[x]]@xpos)>0))
+			hlpEff <- hlp[!(hlp%in%hlpAnno)]
+			hlpEffIndex <- unlist(lapply(hlpEff,function(x) length(overlayInfo[[x]]@start)>0))
+			hlpIndex <- c(hlpAnnoIndex,hlpEffIndex)
+			overlayInfoCorr <- lapply(hlp,function(x) overlayInfo[[x]])
+			gdPlot(trackInfo,minBase=minBase,maxBase=maxBase,overlay=overlayInfoCorr[hlpIndex])
+		}
+	}
+#	if (getWfmFitMethod(object)=="compareGroupsTime" | getWfmFitMethod(object)=="compareGroupsFactor")
+	else if (inherits(inf,"WfmInfCompare"))
+	{
+		if (length(tracks)==1 & tracks[1]==0)
+		{
+			
+			effectsToPlot <- rep(0,noGroups)
+			effectId <- c(1,3,6,10,15)
+			effectNames <- c("2-1","3-2","4-3","5-4","6-5")
+			effectNo <- 1
+			while (effectNo<length(effectsToPlot))
+			{
+				trackInfo[[trackCount]] <- makeGenericArray(intensity=as.matrix(effects[effectId[effectNo]+1,sta:end]),probeStart=Gloc[sta:end],dp=DisplayPars(color="black",ylim=c(min(-1.3,range(effects[effectId[effectNo]+1,sta:end])[1]),max(1.3,range(effects[effectId[effectNo]+1,sta:end])[2]))+c(-0.2,0.2),pointSize=.3,pch=1,lwd=1,type="line"))
+				names(trackInfo)[trackCount] <- effectNames[effectNo]
+				overlayInfo[[trackCount]] <- makeNewTranscriptRectangleOverlay(sigRegions=as.matrix(data.frame(start(regions[[effectId[effectNo]+1]]),end(regions[[effectId[effectNo]+1]]))),location=Gloc,start=sta,end=end,region=c(trackCount,trackCount),dp=DisplayPars(color="darkgrey",alpha=.1))
+				effectNo <- effectNo + 1
+				trackCount <- trackCount + 1
+			}
+			trackInfo[[trackCount]] <- makeGenericArray(intensity=as.matrix(effects[effectId[length(effectsToPlot)-2]+2,sta:end]),probeStart=Gloc[sta:end],dp=DisplayPars(color="black",ylim=c(min(-1.3,range(effects[effectId[length(effectsToPlot)-2]+2,sta:end])[1]),max(1.3,range(effects[effectId[length(effectsToPlot)-2]+2,sta:end])[2]))+c(-0.2,0.2),pointSize=.3,pch=1,lwd=1,type="line"))
+			names(trackInfo)[trackCount] <- "Last-First"
+			overlayInfo[[trackCount]] <- makeNewTranscriptRectangleOverlay(sigRegions=as.matrix(data.frame(start(regions[[effectId[length(effectsToPlot)-2]+2]]),end(regions[[effectId[length(effectsToPlot)-2]+2]]))),location=Gloc,start=sta,end=end,region=c(trackCount,trackCount),dp=DisplayPars(color="darkgrey",alpha=.1))
+			hlp <- which(!unlist(lapply(overlayInfo,is.null)))
+			if (two.strand==TRUE)
+			{
+				hlpAnno <- hlp[1:2]
+			} else
+			{
+				hlpAnno <- hlp[1]
+			}
+			hlpAnnoIndex <- unlist(lapply(hlpAnno,function(x) length(overlayInfo[[x]]@xpos)>0))
+			hlpEff <- hlp[!(hlp%in%hlpAnno)]
+			hlpEffIndex <- unlist(lapply(hlpEff,function(x) length(overlayInfo[[x]]@start)>0))
+			hlpIndex <- c(hlpAnnoIndex,hlpEffIndex)
+			overlayInfoCorr <- lapply(hlp,function(x) overlayInfo[[x]])
+			gdPlot(trackInfo,minBase=minBase,maxBase=maxBase,overlay=overlayInfoCorr[hlpIndex])
+			
+		}
+		if (length(tracks)>1 | tracks[1]!=0)
+		{
+			if (length(tracks)>8)
+			{
+				stop("Too many tracks specified")
+			}
+			groupNames <- sapply(1:noGroups,function(x) paste("Gr",x,sep=""))
+			#groupNames <- sapply(8:13,function(x) paste("D",x,sep=""))
+			firstId <- rep(2:noGroups,1:(noGroups-1))
+			lastId <- sequence(1:(noGroups-1))
+			effectNames <- sapply(1:(noGroups*(noGroups-1)/2),function(x) paste(groupNames[firstId[x]],"-",groupNames[lastId[x]],sep=""))
+			for (i in tracks)
+			{
+				#trackInfo[[trackCount]] <- makeGenericArray(intensity=as.matrix(effects[i+1,sta:end]),probeStart=Gloc[sta:end],dp=DisplayPars(color="black",ylim=c(min(-1.3,range(effects[i+1,sta:end])[1]),max(1.3,range(effects[i+1,sta:end])[2]))+c(-0.2,0.2),pointSize=.3,pch=1,lwd=1,type="line"))
+				trackInfo[[trackCount]] <- makeGenericArray(intensity=as.matrix(effects[i+1,sta:end]),probeStart=Gloc[sta:end],dp=DisplayPars(color="black",ylim=range(effects[,sta:end])+c(-0.2,0.2),pointSize=.3,pch=1,lwd=1,type="line"))
+				names(trackInfo)[trackCount] <- effectNames[i]
+				overlayInfo[[trackCount]] <- makeNewTranscriptRectangleOverlay(sigRegions=as.matrix(data.frame(start(regions[[i+1]]),end(regions[[i+1]]))),location=Gloc,start=sta,end=end,region=c(trackCount,trackCount),dp=DisplayPars(color="darkgrey",alpha=.1))
+				trackCount <- trackCount + 1
+			}
+			hlp <- which(!unlist(lapply(overlayInfo,is.null)))
+			if (two.strand==TRUE)
+			{
+				hlpAnno <- hlp[1:2]
+			} else
+			{
+				hlpAnno <- hlp[1]
+			}
+			hlpAnnoIndex <- unlist(lapply(hlpAnno,function(x) length(overlayInfo[[x]]@xpos)>0))
+			hlpEff <- hlp[!(hlp%in%hlpAnno)]
+			hlpEffIndex <- unlist(lapply(hlpEff,function(x) length(overlayInfo[[x]]@start)>0))
+			hlpIndex <- c(hlpAnnoIndex,hlpEffIndex)
+			overlayInfoCorr <- lapply(hlp,function(x) overlayInfo[[x]])
+			gdPlot(trackInfo,minBase=minBase,maxBase=maxBase,overlay=overlayInfoCorr[hlpIndex])
+		}
+		
+	}
+#	if (getWfmFitMethod(object)=="meansByGroupTime" | getWfmFitMethod(object)=="meansByGroupFactor")
+	else if (inherits(inf,"WfmInfMeans"))
+	{
+		if (length(tracks)==1 & tracks[1]==0)
+		{
+			effectNames <- sapply(1:noGroups,function(x) paste("Gr",x,sep=""))
+			for (i in 1:noGroups)
+			{
+				trackInfo[[trackCount]] <- makeGenericArray(intensity=as.matrix(effects[i,sta:end]),probeStart=Gloc[sta:end],dp=DisplayPars(color="black",ylim=c(min(-1.3,range(effects[i,sta:end])[1]),max(1.3,range(effects[i,sta:end])[2]))+c(-0.2,0.2),pointSize=.3,pch=1,lwd=1,type="line"))
+				names(trackInfo)[trackCount] <- effectNames[i]
+				overlayInfo[[trackCount]] <- makeNewTranscriptRectangleOverlay(sigRegions=as.matrix(data.frame(start(regions[[i]]),end(regions[[i]]))),location=Gloc,start=sta,end=end,region=c(trackCount,trackCount),dp=DisplayPars(color="darkgrey",alpha=.1))
+				trackCount <- trackCount + 1
+			}
+			hlp <- which(!unlist(lapply(overlayInfo,is.null)))
+			if (two.strand==TRUE)
+			{
+				hlpAnno <- hlp[1:2]
+			} else
+			{
+				hlpAnno <- hlp[1]
+			}
+			hlpAnnoIndex <- unlist(lapply(hlpAnno,function(x) length(overlayInfo[[x]]@xpos)>0))
+			hlpEff <- hlp[!(hlp%in%hlpAnno)]
+			hlpEffIndex <- unlist(lapply(hlpEff,function(x) length(overlayInfo[[x]]@start)>0))
+			hlpIndex <- c(hlpAnnoIndex,hlpEffIndex)
+			overlayInfoCorr <- lapply(hlp,function(x) overlayInfo[[x]])
+			gdPlot(trackInfo,minBase=minBase,maxBase=maxBase,overlay=overlayInfoCorr[hlpIndex])
+		}
+		if (length(tracks)>1 | tracks[1]!=0)
+		{
+			if (length(tracks)>8)
+			{
+				stop("Too many tracks specified")
+			}
+			#effectNames <- sapply(1:noGroups,function(x) paste("Gr",x,sep=""))
+			effectNames <- paste("day",8:13,sep="")
+			for (i in tracks)
+			{
+				trackInfo[[trackCount]] <- makeGenericArray(intensity=as.matrix(effects[i,sta:end]),probeStart=Gloc[sta:end],dp=DisplayPars(color="black",ylim=c(min(-1.3,range(effects[i,sta:end])[1]),max(1.3,range(effects[i,sta:end])[2]))+c(-0.2,0.2),pointSize=.3,pch=1,lwd=1,type="line"))
+				#trackInfo[[trackCount]] <- makeGenericArray(intensity=as.matrix(effects[i,sta:end]),probeStart=Gloc[sta:end],dp=DisplayPars(color="black",ylim=range(effects[,sta:end])+c(-0.2,0.2),pointSize=.3,pch=1,lwd=1,type="line"))
+				names(trackInfo)[trackCount] <- effectNames[i]
+				overlayInfo[[trackCount]] <- makeNewTranscriptRectangleOverlay(sigRegions=as.matrix(data.frame(start(regions[[i]]),end(regions[[i]]))),location=Gloc,start=sta,end=end,region=c(trackCount,trackCount),dp=DisplayPars(color="darkgrey",alpha=.1))
+				trackCount <- trackCount + 1
+			}
+			hlp <- which(!unlist(lapply(overlayInfo,is.null)))
+			if (two.strand==TRUE)
+			{
+				hlpAnno <- hlp[1:2]
+			} else
+			{
+				hlpAnno <- hlp[1]
+			}
+			hlpAnnoIndex <- unlist(lapply(hlpAnno,function(x) length(overlayInfo[[x]]@xpos)>0))
+			hlpEff <- hlp[!(hlp%in%hlpAnno)]
+			hlpEffIndex <- unlist(lapply(hlpEff,function(x) length(overlayInfo[[x]]@start)>0))
+			hlpIndex <- c(hlpAnnoIndex,hlpEffIndex)
+			overlayInfoCorr <- lapply(hlp,function(x) overlayInfo[[x]])
+			gdPlot(trackInfo,minBase=minBase,maxBase=maxBase,overlay=overlayInfoCorr[hlpIndex])
+		}
+	}
+#	if (getWfmFitMethod(object)=="effectsTime")
+	else if (inherits(fit,"WfmFitTime") & inherits(inf,"WfmInfEffects"))
+	{
+		fct <- 0.2390457
+		if (tracks==0 | (length(tracks)==2 & tracks[1]==1 & tracks[2]==2) | tracks==1 | tracks==2)
+		{
+			if (tracks==0)
+			{
+				tracks <- c(1,2)
+			}
+			effectNames <- c("Linear","Quadratic")
+			for (i in tracks)
+			{
+				#trackInfo[[trackCount]] <- makeGenericArray(intensity=as.matrix(effects[i+1,sta:end]),probeStart=Gloc[sta:end],dp=DisplayPars(color="black",ylim=range(effects[i+1,sta:end])+c(-0.5,0.5),pointSize=.3,pch=1,lwd=1,type="line"))
+				trackInfo[[trackCount]] <- makeGenericArray(intensity=as.matrix(effects[i+1,sta:end]*fct),probeStart=Gloc[sta:end],dp=DisplayPars(color="black",ylim=range(effects[i+1,sta:end]*fct)+c(-0.2,0.2),pointSize=.3,pch=1,lwd=1,type="line"))
+				names(trackInfo)[trackCount] <- effectNames[i]
+				overlayInfo[[trackCount]] <- makeNewTranscriptRectangleOverlay(sigRegions=as.matrix(data.frame(start(regions[[i+1]]),end(regions[[i+1]]))),location=Gloc,start=sta,end=end,region=c(trackCount,trackCount),dp=DisplayPars(color="darkgrey",alpha=.1))
+				trackCount <- trackCount + 1
+			}
+			hlp <- which(!unlist(lapply(overlayInfo,is.null)))
+			if (two.strand==TRUE)
+			{
+				hlpAnno <- hlp[1:2]
+			} else
+			{
+				hlpAnno <- hlp[1]
+			}
+			hlpAnnoIndex <- unlist(lapply(hlpAnno,function(x) length(overlayInfo[[x]]@xpos)>0))
+			hlpEff <- hlp[!(hlp%in%hlpAnno)]
+			hlpEffIndex <- unlist(lapply(hlpEff,function(x) length(overlayInfo[[x]]@start)>0))
+			hlpIndex <- c(hlpAnnoIndex,hlpEffIndex)
+			overlayInfoCorr <- lapply(hlp,function(x) overlayInfo[[x]])
+			gdPlot(trackInfo,minBase=minBase,maxBase=maxBase,overlay=overlayInfoCorr[hlpIndex])
+		} else
+		{
+				stop("tracks argument has been misspecified")
+		}
+	}
+	else if (inherits(inf,"WfmInfCustom")) {
+		stop("Custom Inference plot not implemented!")
+	}
+	else {
+		stop("Unknown Design")	  
+	}
+}
\ No newline at end of file

Modified: pkg/R/initialize-methods.R
===================================================================
--- pkg/R/initialize-methods.R	2011-12-21 09:38:39 UTC (rev 17)
+++ pkg/R/initialize-methods.R	2012-02-15 16:55:39 UTC (rev 18)
@@ -9,7 +9,7 @@
 })
 
 #WfmFit
-setMethod("initialize","Wfm",function(.Object,betaMAP,varbetaMAP,smoothPar,varEps,dataOrigSpace,dataWaveletSpace,design,phenoData,method,genome.info,n.levels,probePosition,wave.filt,Kj,prior,alpha,delta,two.sided,rescale,sigProbes,regions,GlocRegions,FDR,CI,F,varF,eff,varEff)
+setMethod("initialize","WfmFit",function(.Object,betaMAP,varbetaMAP,smoothPar,varEps,dataOrigSpace,dataWaveletSpace,design.matrix,phenoData,genome.info,n.levels,probePosition,wave.filt,Kj,prior,F,varF,P,Z,noGroups,replics)
 {
 	.Object at betaMAP <- betaMAP
 	.Object at varbetaMAP <- varbetaMAP
@@ -17,32 +17,42 @@
 	.Object at varEps <- varEps
 	.Object at dataOrigSpace <- dataOrigSpace
 	.Object at dataWaveletSpace <- dataWaveletSpace
-	.Object at design <- design
+	.Object at design.matrix <- design.matrix
 	.Object at phenoData <- phenoData
-	.Object at method <- method
+#	.Object at design <- design
 	.Object at genome.info <- genome.info
 	.Object at n.levels <- n.levels
 	.Object at probePosition <- probePosition
 	.Object at wave.filt <- wave.filt
 	.Object at Kj <- Kj
 	.Object at prior <- prior
+#	.Object at rescale <- rescale
+	.Object at F <- F
+	.Object at varF <- varF
+	.Object at P <- P
+	.Object at Z <- Z
+	.Object at noGroups <- noGroups
+	.Object at replics <- replics
+	return(.Object)
+})
+
+
+# WfmInf
+setMethod("initialize","WfmInf",function(.Object, alpha, delta, two.sided, sigProbes, regions, GlocRegions, FDR, CI, eff, varEff)
+{
 	.Object at alpha <- alpha
 	.Object at delta <- delta
 	.Object at two.sided <- two.sided
-	.Object at rescale <- rescale
 	.Object at sigProbes <- sigProbes
 	.Object at regions <- regions
 	.Object at GlocRegions <- GlocRegions
 	.Object at FDR <- FDR
 	.Object at CI <- CI
-	.Object at F <- F
-	.Object at varF <- varF
-	.Object at eff <- eff
-	.Object at varEff <- varEff
+ 	.Object at eff <- eff
+ 	.Object at varEff <- varEff
 	return(.Object)
 })
 
-
 #genomeInfo
 setMethod("initialize","genomeInfo",function(.Object,chromosome,strand,minPos,maxPos)
 {
@@ -54,9 +64,63 @@
 })
 
 #waveTilingFeatureSet
-setMethod("initialize","WaveTilingFeatureSet",function (.Object, ...)
+setMethod("initialize","WaveTilingFeatureSet",function (.Object)
 {
 	callNextMethod(.Object);
 }
 )
 
+## Subclasses
+setMethod("initialize","WfmFitTime",function (.Object, ...)
+{
+	callNextMethod(.Object, ...);
+}
+)
+
+setMethod("initialize","WfmFitCircadian",function (.Object, ...)
+{
+	callNextMethod(.Object, ...);
+}
+)
+
+setMethod("initialize","WfmFitFactor",function (.Object, ...)
+{
+	callNextMethod(.Object, ...);
+}
+)
+
+setMethod("initialize","WfmFitCustom",function (.Object, ...)
+{
+	callNextMethod(.Object, ...);
+}
+)
+
+setMethod("initialize","WfmInfCustom",function (.Object, ...)
+{
+	callNextMethod(.Object, ...);
+}
+)
+
+setMethod("initialize","WfmInfCompare",function (.Object, ...)
+{
+	callNextMethod(.Object, ...);
+}
+)
+
+setMethod("initialize","WfmInfEffects",function (.Object, ...)
+{
+	callNextMethod(.Object, ...);
+}
+)
+
+setMethod("initialize","WfmInfMeans",function (.Object, ...)
+{
+	callNextMethod(.Object, ...);
+}
+)
+
+setMethod("initialize","WfmInfOverallMean",function (.Object, ...)
+{
+	callNextMethod(.Object, ...);
+}
+)

Modified: pkg/R/methods-WaveTilingFeatureSet.R
===================================================================
--- pkg/R/methods-WaveTilingFeatureSet.R	2011-12-21 09:38:39 UTC (rev 17)
+++ pkg/R/methods-WaveTilingFeatureSet.R	2012-02-15 16:55:39 UTC (rev 18)
@@ -79,7 +79,7 @@
 
 setMethod("filterOverlap",signature("WaveTilingFeatureSet"),function(object,remap=TRUE,fastaFile,chrId,strand=c("forward","reverse","both"),MM=FALSE)
 {
-	if (class(object)!="WaveTilingFeatureSet")
+	if (!inherits(object,"WaveTilingFeatureSet")) #class(object)!="WaveTilingFeatureSet")
 	{
 		stop("class of object is not WaveTilingFeatureSet.")
 	}
@@ -254,7 +254,7 @@
 
 setMethod("selectProbesFromTilingFeatureSet",signature("WaveTilingFeatureSet"),function(object,chromosome,strand=c("forward","reverse"),minPos,maxPos)
 {
-	if (class(object)!="TilingFeatureSet")
+	if (!inherits(object,"TilingFeatureSet")) #class(object)!="TilingFeatureSet")
 		{
 			stop("class of object is not TilingFeatureSet.")
 		}
@@ -308,86 +308,17 @@
 })
 
 
-setMethod("makeDesign",signature("WaveTilingFeatureSet"),function(object,method=c("twoGroup","compareGroupsTime","compareGroupsFactor","circadian","twoFactors","meansByGroupTime","meansByGroupFactor","effectsTime"),factor.levels=NULL)
-{
-	replics <- getReplics(object)
-	noGroups <- getNoGroups(object)
-	if (method=="twoGroup")
-	{
-		Xorig <- matrix(0,nrow=dim(pData(object))[1],ncol=2)
-		Xorig[,1] <- 1
-		Xorig[,2] <- rep(c(1,-1),replics)
-	}
-	if (method=="compareGroupsTime" | method=="meansByGroupTime" | method=="effectsTime")
-	{
-		Xorig <- matrix(0,nrow=dim(pData(object))[1],ncol=noGroups)
-		orderedFactor <- factor(1:noGroups,ordered=TRUE)
-		desPoly <- lm(rnorm(noGroups)~orderedFactor,x=TRUE)$x
-		Xorig[,1] <- 1
-		Xorig[,2:noGroups] <- apply(desPoly[,2:noGroups],2,rep,replics)
-		
-	}
-	if (method=="circadian")
-	{
-		Xorig <- matrix(0,nrow=dim(pData(object))[1],ncol=3)
-		Xorig[,1] <- 1
-		Xorig[,2] <- rep(sin(seq(0,2*pi-(pi/noGroups),2*pi/noGroups)),replics)
-		Xorig[,3] <- rep(cos(seq(0,2*pi-(pi/noGroups),2*pi/noGroups)),replics)
-	}
-	if (method=="compareGroupsFactor" | method=="meansByGroupFactor")
-	{
-		Xorig <- matrix(0,nrow=dim(pData(object))[1],ncol=noGroups)
-		desHelmert <- contr.helmert(noGroups)
-		Xorig[,1] <- 1
-		Xorig[,2:noGroups] <- apply(desHelmert[,1:(noGroups-1)],2,rep,replics)
-	}
-	if (method=="twoFactors")
-	{
-		if (is.null(factor.levels))
-		{
-			stop("No proper factor levels given.")
-		}
-		Xorig <- matrix(0,nrow=sum(replics),ncol=prod(factor.levels))
-		Xorig[,1] <- 1
-		desTrt1 <- contr.treatment(factor.levels[1])
-		desTrt2 <- contr.treatment(factor.levels[2])
-		#desHelmert1 <- contr.helmert(factor.levels[1])
-		#desHelmert2 <- contr.helmert(factor.levels[2])
-		replicsMat <- matrix(replics,nrow=factor.levels[1],ncol=factor.levels[2],byrow=TRUE)
-		desFac1 <- apply(desTrt1,2,function(x) unlist(sapply(1:nrow(replicsMat),function(y) rep(x[y],sum(replicsMat[y,])))))
-		Xorig[,2:(factor.levels[1])] <- desFac1
-		desFac2 <- apply(desTrt2,2,function(x) unlist(sapply(1:nrow(replicsMat),function(y) rep(x,replicsMat[y,]))))
-		Xorig[,(factor.levels[1]+1):(factor.levels[1]+factor.levels[2]-1)] <- desFac2
-		# include interactions
-		colno <- sum(factor.levels)
-		for (i in 1:ncol(desFac1))
-		{
-			for (j in 1:ncol(desFac2))
-			{
-				Xorig[,colno] <- desFac1[,i]*desFac2[,j]
-				colno <- colno + 1
-			}
-		}
-	}
-	out <- list()
-	X.qr <- qr(Xorig)
-	out$Xorthnorm <- qr.Q(X.qr)
-	out$Xorig <- Xorig
-	return(out)
-})
 
-
[TRUNCATED]

To get the complete diff run:
    svnlook diff /svnroot/wavetiling -r 18


More information about the Wavetiling-commits mailing list