R/reformatmenu.R

Defines functions peaksDetector unionReads myfindOverlaps

#oneChannelGUI reformatGdl converts a gdl object in an expression matrix containing NGS counts
#              internal functions:
#                           .unionReads
#                           .peaksDetector
#                           .myfindOverlaps
#oneChannelGUI rankingConversion converts intensity signal in normalized ranks
#               in case of NGS counts ranking, since ranking goeas from 0 (low counts) to 1 (high counts) is  moltipled for 
#               16 to give a range of signal between 0 and 16
#This function is in normalized menu
#oneChannelGUI log2Conversion makes a log2 transformation of the counts
################################################################################
"reformatGdl" <- function(){
          Try(tkconfigure(.affylmGUIglobals$ttMain,cursor="watch"))
   #       require(chipseq) || stop("\nMissing chipseq library\n")
          #error if no data are loaded
          Try(NGS.Available <- get("NGS.Available",envir=affylmGUIenvironment))
          if(!NGS.Available){
              Try(tkmessageBox(title="Reformatting NGS",message="No raw NGS data have been loaded.	Please try New or Open from the File menu.",type="ok",icon="error"))
				    	Try(tkfocus(.affylmGUIglobals$ttMain))
              Try(return())
          }
          ###########################################
          Try(ttGetExtension<-tktoplevel(.affylmGUIglobals$ttMain))
          Try(tkwm.deiconify(ttGetExtension))
          Try(tkgrab.set(ttGetExtension))
          Try(tkfocus(ttGetExtension))
          Try(tkwm.title(ttGetExtension,"Extending reads length"))
          Try(tkgrid(tklabel(ttGetExtension,text="    ")))
          Try(ExtensionText <- "35")
          Try(Local.Extension <- tclVar(init=ExtensionText))
          Try(entry.Extension <-tkentry(ttGetExtension,width="10",font=.affylmGUIglobals$affylmGUIfont2,textvariable=Local.Extension,bg="white"))
          Try(tkgrid(tklabel(ttGetExtension,text="Please enter the size of read extension.",font=.affylmGUIglobals$affylmGUIfont2)))
          Try(tkgrid(entry.Extension))
          onOK <- function()
          {
              Try(ExtensionText <- tclvalue(Local.Extension))
              if(nchar(ExtensionText)==0)
                ExtensionText <- "Unselected"
                Try(assign("ExtensionText",ExtensionText,affylmGUIenvironment))
                Try(tclvalue(.affylmGUIglobals$ExtensionTcl) <- ExtensionText)
                Try(tkgrab.release(ttGetExtension));Try(tkdestroy(ttGetExtension));Try(tkfocus(.affylmGUIglobals$ttMain))
          }
          Try(OK.but <-tkbutton(ttGetExtension,text="   OK   ",command=onOK,font=.affylmGUIglobals$affylmGUIfont2))
          Try(tkgrid(tklabel(ttGetExtension,text="    ")))
          Try(tkgrid(OK.but))
          Try(tkgrid.configure(OK.but))
          Try(tkgrid(tklabel(ttGetExtension,text="       ")))
          Try(tkfocus(entry.Extension))
          Try(tkbind(entry.Extension, "<Return>",onOK))
          Try(tkbind(ttGetExtension, "<Destroy>", function(){Try(tkgrab.release(ttGetExtension));Try(tkfocus(.affylmGUIglobals$ttMain));return(0)}))
          Try(tkwait.window(ttGetExtension))
          Try(tkfocus(.affylmGUIglobals$ttMain))
          Try(extensionSize <- as.numeric(get("ExtensionText",envir=affylmGUIenvironment)))
          #extending reads length
          #Try(extRanges <- gdapply(gdl, extendReads, seqLen = extensionSize))
          ###########################################
          #defining the threshold for peak definition
          Try(ttGetreadThr<-tktoplevel(.affylmGUIglobals$ttMain))
          Try(tkwm.deiconify(ttGetreadThr))
          Try(tkgrab.set(ttGetreadThr))
          Try(tkfocus(ttGetreadThr))
          Try(tkwm.title(ttGetreadThr,"Minimal number of read in a peak"))
          Try(tkgrid(tklabel(ttGetreadThr,text="    ")))
          Try(readThrText <- "8")
          Try(Local.readThr <- tclVar(init=readThrText))
          Try(entry.readThr <-tkentry(ttGetreadThr,width="4",font=.affylmGUIglobals$affylmGUIfont2,textvariable=Local.readThr,bg="white"))
          Try(tkgrid(tklabel(ttGetreadThr,text="Please enter the minimal number of reads needed \nto define a peak.",font=.affylmGUIglobals$affylmGUIfont2)))
          Try(tkgrid(entry.readThr))
          onOK <- function()
          {
              Try(readThrText <- tclvalue(Local.readThr))
              if(nchar(readThrText)==0)
                readThrText <- "Unselected"
                Try(assign("readThrText",readThrText,affylmGUIenvironment))
                Try(tclvalue(.affylmGUIglobals$readThrTcl) <- readThrText)
                Try(tkgrab.release(ttGetreadThr));Try(tkdestroy(ttGetreadThr));Try(tkfocus(.affylmGUIglobals$ttMain))
          }
          Try(OK.but <-tkbutton(ttGetreadThr,text="   OK   ",command=onOK,font=.affylmGUIglobals$affylmGUIfont2))
          Try(tkgrid(tklabel(ttGetreadThr,text="    ")))
          Try(tkgrid(OK.but))
          Try(tkgrid.configure(OK.but))
          Try(tkgrid(tklabel(ttGetreadThr,text="       ")))
          Try(tkfocus(entry.readThr))
          Try(tkbind(entry.readThr, "<Return>",onOK))
          Try(tkbind(ttGetreadThr, "<Destroy>", function(){Try(tkgrab.release(ttGetreadThr));Try(tkfocus(.affylmGUIglobals$ttMain));return(0)}))
          Try(tkwait.window(ttGetreadThr))
          Try(tkfocus(.affylmGUIglobals$ttMain))
          Try(readsMinSize <- as.numeric(get("readThrText",envir=affylmGUIenvironment)))
          #union at chr level of all reads
          Try(gdl <- get("NGSdata",envir=affylmGUIenvironment))
          Try(readsMinSize <- readsMinSize * length(gdl))
          Try(chrs.length <- get("chrsLength",envir=affylmGUIenvironment))
          #now I have to start to count reads from any experiment in each peak
          Try(mymat <- .unionReads(gdl, chrs.length, seqLen=extensionSize, lower=readsMinSize))
          ######################################################################
          #Creating the expression set
          Try(Targets <- get("Targets",envir=affylmGUIenvironment))
          Try(rownames(Targets) <- Targets$FileName) 
          Try(var.tmp.pd<-data.frame(names(Targets)))
          Try(names(var.tmp.pd)<-"labelDescription" )
          Try(rownames(var.tmp.pd)<-names(Targets))
          Try(tmp1.pd<-new("AnnotatedDataFrame", data=Targets, varMetadata=var.tmp.pd))
          Try(NormalizedAffyData <- new("ExpressionSet", exprs= mymat , phenoData=tmp1.pd, annotation=""))
	        Try(tkconfigure(.affylmGUIglobals$ttMain,cursor="arrow"))
          Try(assign("NormalizedAffyData.Available",TRUE,affylmGUIenvironment))
	        Try(assign("NormalizedAffyData",NormalizedAffyData,affylmGUIenvironment))
          Try(tkdelete(.affylmGUIglobals$mainTree,"NormalizedAffyData.Status"))
	        Try(tkinsert(.affylmGUIglobals$mainTree,"end","NormalizedAffyData","NormalizedAffyData.Status" ,text="NGS counts",font=.affylmGUIglobals$affylmGUIfontTree))
}
################################################################################
# used in reformatGdi by .unionReads
.peaksDetector <- function(x, seqLen, strand, width, lower){
                          g <- extendReads(x, seqLen = seqLen, strand = strand)
                          g1 <- coverage(g, width = width)
                          s <- slice(g1, lower = lower)
                          peaks.ir <- IRanges(start = start(s), end = end(s))
                          return(peaks.ir)
}

################################################################################
#used by reformatGdi
#gdl a GenomeDataList object
#chrs.length size of each chr
#seqLen size of extesion of reads
#lower numer of reads got by chance 
.unionReads <- function(gdl, chrs.length, seqLen , lower){
     #     Try(require(BSgenome) || stop("\nPackage BSgenome is missing.\n"))
          Try(union.list <- list())
          for (j in 1:length(gdl[[1]])){
             Try(n.reads <- NULL)
             Try(p.reads <- NULL)
             for (i in 1:length(gdl)){
                 Try(n.reads <- c(n.reads, gdl[[i]][[j]][["-"]]))
                 Try(p.reads <- c(p.reads, gdl[[i]][[j]][["+"]])) 
             }
             Try(union.list[[names(gdl[[1]][j])]] <-  list("-" = n.reads, "+" = p.reads))
          }
          Try(union.gd <- GenomeData(listData = union.list, organism = organism(gdl[[1]]), provider = provider(gdl[[1]]), providerVersion = providerVersion(gdl[[1]])))
          Try(union.list <- list())
          Try(union.p <- list())
          Try(union.n <- list())
          
          for(i in 1:length(union.gd)){
                Try(tmp.p <- gdapply(union.gd[i], .peaksDetector, seqLen = seqLen, strand = "+", width = chrs.length[i], lower = lower))
                Try(tmp1.p <- cbind(start(tmp.p[[1]]),end(tmp.p[[1]])))
                Try(union.p[[names(union.gd[i])]] <- tmp1.p)
                Try(tmp.n <- gdapply(union.gd[i], .peaksDetector, seqLen = seqLen, strand = "-", width = chrs.length[i], lower = lower))
                Try(tmp1.n <- cbind(start(tmp.n[[1]]),end(tmp.n[[1]])))
                Try(union.n[[names(union.gd[i])]] <- tmp1.n)

          }
          Try(union.list <- list("-" = union.n, "+" = union.p))
          Try(mymat <- NULL)
          Try(cat("\nStart reformatting\n"))
          #loop on all chrs c is chr
          for(c in 1:length(union.list[["-"]])){
             # loop on all sample s is sample
             Try(mydf.m <- NULL)
             Try(mydf.p <- NULL) 
             for(s in 1:length(gdl)){
               Try(cat("reformatting sample", names(gdl[s]), " ", names(gdl[[s]][c]), "\n"))
               if(identical(names(union.list[["-"]][c]), names(gdl[[s]][c]))){
                 Try(myquery.m <- IRanges(start = union.list[["-"]][[c]][,1], end = union.list[["-"]][[c]][,2]))
                 Try(mysubject.m <- IRanges(start = gdl[[s]][[c]][["-"]], width = rep(1, length(gdl[[s]][[c]][["-"]]))))
                 Try(peaks.m <- sapply(myquery.m, .myfindOverlaps , y=mysubject.m))
                 Try(peaks.m.length <- sapply(peaks.m, length))
                 Try(mydf.m <- cbind(mydf.m, peaks.m.length))

                Try(myquery.p <- IRanges(start = union.list[["+"]][[c]][,1], end = union.list[["+"]][[c]][,2]))
                Try(mysubject.p <- IRanges(start = gdl[[s]][[c]][["+"]], width = rep(1, length(gdl[[s]][[c]][["+"]]))))
                Try(peaks.p <- sapply(myquery.p, .myfindOverlaps, y=mysubject.p))
                Try(peaks.p.length <- sapply(peaks.p, length))
                Try(mydf.p <- cbind(mydf.p, peaks.p.length))
              } else {
                  Try(cat("\nThere are some problems in the correct association between chrs during the reformatting of NGS data in an expression Matrix."))
                  Try(cat("\nPlease contact oneChannelGUI mantainer indicating the error NGS500\n"))
                  Try(return(NULL))
              }
            }
            if(dim(mydf.m)[1] > 0 || dim(mydf.p)[1] > 0){
               Try(dimnames(mydf.m) <- list(paste(names(union.list[["-"]][c]), "-", union.list[["-"]][[c]][,1], union.list[["-"]][[c]][,2], sep="_"), names(gdl)))
               Try(dimnames(mydf.p) <- list(paste(names(union.list[["+"]][c]), "+", union.list[["+"]][[c]][,1], union.list[["+"]][[c]][,2], sep="_"), names(gdl)))
               Try(mymat <- rbind(mymat, mydf.m, mydf.p))
            }
         }
         Try(cat("\nEnd reformatting\n"))
         Try(return(mymat))
}
################################################################################
.myfindOverlaps <- function(x,y) findOverlaps(x, y)
#################################################################################
"log2Conversion" <- function(){
  if(affylmGUIenvironment$NormalizedAffyData.Available){
                     Try(NormalizedAffyData <- get("NormalizedAffyData", env=affylmGUIenvironment))
  } else Try(tkmessageBox(title="NGS data set normalization",message="Normalized data are not available.", type="ok", icon="error"))
  Try(tkconfigure(.affylmGUIglobals$ttMain,cursor="watch"))
  Try(tmp.rank <- esApply(NormalizedAffyData, 2, log2))
  if(affylmGUIenvironment$NGS.Available){
        Try(exprs(NormalizedAffyData) <- tmp.rank)
  } else Try(exprs(NormalizedAffyData) <- tmp.rank)
  Try(assign("NormalizedAffyData", NormalizedAffyData,  env=affylmGUIenvironment))
  Try(assign("NGSconversion.available", TRUE,  env=affylmGUIenvironment))
  tkdelete(.affylmGUIglobals$mainTree,"NormalizedAffyData.Status")
  tkinsert(.affylmGUIglobals$mainTree,"end","NormalizedAffyData","NormalizedAffyData.Status" ,text="Available log2 transformed counts",font=.affylmGUIglobals$affylmGUIfontTree)
  Try(tkmessageBox(title="Large data set normalization",message="Data are now available as log2 transformed counts in NormalizedAffyData.", type="ok")) 
Try(tkconfigure(.affylmGUIglobals$ttMain,cursor="arrow"))
}

Try the oneChannelGUI package in your browser

Any scripts or data that you put into this service are public.

oneChannelGUI documentation built on Nov. 17, 2017, 11:02 a.m.