aboutsummaryrefslogtreecommitdiff
path: root/scripts/ctl_analysis.R
diff options
context:
space:
mode:
authorAlexander Kabui2022-01-14 02:42:03 +0300
committerBonfaceKilz2022-01-22 09:23:14 +0300
commit3ab0c245d2bbfe695b5892e5a1a11c78b2cbce50 (patch)
treea8fc07c0d944c26300ed69469e0ffffdba51f627 /scripts/ctl_analysis.R
parent827811810bd37aaca7456b74efc80779f63dcc04 (diff)
downloadgenenetwork3-3ab0c245d2bbfe695b5892e5a1a11c78b2cbce50.tar.gz
function override:fix target specific output file for network
Diffstat (limited to 'scripts/ctl_analysis.R')
-rw-r--r--scripts/ctl_analysis.R116
1 files changed, 114 insertions, 2 deletions
diff --git a/scripts/ctl_analysis.R b/scripts/ctl_analysis.R
index 4002adb..ff7aa94 100644
--- a/scripts/ctl_analysis.R
+++ b/scripts/ctl_analysis.R
@@ -97,12 +97,124 @@ for (trait in phenoData$trait_db_list)
# rename coz of duplicate key names
-network_file_path = file.path(input$imgDir,paste("ctlnet","random",".sif",sep=""))
+network_file_path = file.path(input$imgDir,paste("ctlnet","randomz",sep=""))
file.create(network_file_path)
-ctl_network = CTLnetwork(ctls, significance = 0.05, LODdrop = 2,short = FALSE, add.qtls = FALSE, file = network_file_path, verbose = TRUE)
+
+
+
+# temp fix override ctlnetwork function to target certain output file
+
+
+CTLnetworkGn<- function(CTLobject, mapinfo, significance = 0.05, LODdrop = 2, what = c("names","ids"), short = FALSE, add.qtls = FALSE, file = "",tmpDir = "",verbose = TRUE){
+ if(missing(CTLobject) || is.null(CTLobject)) stop("argument 'CTLobject' is missing, with no default")
+ if("CTLscan" %in% class(CTLobject)) CTLobject = list(CTLobject)
+ if(length(what) > 1) what = what[1]
+
+ results <- NULL
+ significant <- CTLsignificant(CTLobject, significance, what = "ids")
+ if(!is.null(significant)){
+ all_m <- NULL; all_p <- NULL;
+ nodefile=""; netfile = "";
+ if(file != ""){
+
+ netfile = file.path(tmpDir,paste("ctlnet",file,".sif",sep=""))
+
+ nodefile = file.path(tmpDir,paste("ctlnet",file,".nodes",sep=""))
+
+ # netfile <- paste("ctlnet",file,".sif",sep="")
+ # nodefile <- paste("ctlnet",file,".nodes",sep="")
+ }
+ cat("",file=netfile); cat("",file=nodefile);
+ if(verbose) cat("NETWORK.SIF\n")
+ edges <- NULL
+ for(x in 1:nrow(significant)){
+ data <- as.numeric(significant[x,])
+ CTLscan <- CTLobject[[data[1]]]
+ markern <- rownames(CTLscan$dcor)
+ traitsn <- colnames(CTLscan$dcor)
+ name <- ctl.name(CTLscan)
+ if(what=="ids"){
+ tid <- which(traitsn %in% ctl.name(CTLobject[[data[1]]]))
+ name <- paste("P",tid,sep="")
+ markern <- paste("M",1:nrow(CTLobject[[data[1]]]$dcor), sep="")
+ traitsn <- paste("P", 1:ncol(CTLobject[[data[1]]]$dcor), sep="")
+ }
+ if(add.qtls){ # Add QTL to the output SIF
+ bfc <- length(CTLscan$qtl)
+ above <- which(CTLscan$qtl > -log10(significance))
+ qtlnms <- names(above); qtlmid <- 1
+ for(m in above){
+ cat(name,"\tQTL\t",markern[m],"\tQTL\t",CTLscan$qtl[m],"\n",sep="",file=netfile,append=TRUE)
+ all_m <- CTLnetwork.addmarker(all_m, mapinfo, markern[data[2]], qtlnms[qtlmid])
+ qtlmid <- qtlmid+1
+ }
+ }
+ lod <- CTLscan$ctl[data[2],data[3]]
+ qlod1 <- CTLscan$qtl[data[2]]
+ qlod2 <- qlod1
+ edgetype <- NA
+ if(length(CTLobject) >= data[3]){ # Edge type based on QTL LOD scores
+ qlod2 <- CTLobject[[data[3]]]$qtl[data[2]]
+ if((qlod1-qlod2) > LODdrop){
+ edgetype <- 1
+ }else if((qlod1-qlod2) < -LODdrop){
+ edgetype <- -1
+ }else{ edgetype <- 0; }
+ } else {
+ cat("Warning: Phenotype", data[3], "from", data[1], "no CTL/QTL information\n")
+ qlod2 <- NA;
+ }
+ #Store the results
+ results <- rbind(results, c(data[1], data[2], data[3], lod, edgetype, qlod1, qlod2))
+
+ if(nodefile == "" && !verbose){ }else{
+ if(short){
+ edge <- paste(name,traitsn[data[3]])
+ edgeI <- paste(traitsn[data[3]],name)
+ if(!edge %in% edges && !edgeI %in% edges){
+ cat(name, "\t", markern[data[2]],"\t", traitsn[data[3]],"\n",file=netfile, append=TRUE,sep="")
+ edges <- c(edges,edge)
+ }
+ }else{
+ cat(name, "\t", "CTL_", data[1],"_",data[3], "\t", markern[data[2]],file=netfile, append=TRUE,sep="")
+ cat("\tCTL\t", lod, "\n", file=netfile, append=TRUE,sep="")
+ cat(markern[data[2]], "\t", "CTL_", data[1],"_",data[3], "\t",file=netfile, append=TRUE,sep="")
+ cat(traitsn[data[3]],"\tCTL\t", lod, "\n", file=netfile,append=TRUE,sep="")
+ }
+ }
+ all_m <- CTLnetwork.addmarker(all_m, mapinfo, markern[data[2]], rownames(CTLscan$dcor)[data[2]])
+ all_p <- unique(c(all_p, name, traitsn[data[3]]))
+ }
+ colnames(results) <- c("TRAIT1","MARKER","TRAIT2","LOD_C","CAUSAL","LOD_T1","LOD_T2")
+ if(verbose) cat("NODE.DESCRIPTION\n")
+ if(nodefile == "" && !verbose){ }else{
+ for(m in all_m){ cat(m,"\n", sep="", file=nodefile, append=TRUE); }
+ for(p in all_p){ cat(p,"\tPHENOTYPE\n", sep="", file=nodefile, append=TRUE); }
+ }
+ }
+ if(!is.null(results)){
+ class(results) <- c(class(results),"CTLnetwork")
+ }
+ invisible(results)
+}
+
+CTLnetwork.addmarker <- function(markers, mapinfo, name, realname){
+ if(!missing(mapinfo)){
+ id <- which(rownames(mapinfo) %in% realname)
+ fname <- paste(name,"\tMARKER\t",mapinfo[id,1],"\t",mapinfo[id,2],sep="")
+ markers <- unique(c(markers, fname))
+ }
+ return(markers)
+}
+
+
+# generating network
+
+
+ctl_network = CTLnetworkGn(ctls, significance = 0.05, LODdrop = 2,short = FALSE, add.qtls = FALSE,tmpDir=input$imgDir,file = "random", verbose = TRUE)