659 lines
36 KiB
R
659 lines
36 KiB
R
#!/usr/bin/env Rscript
|
|
# This R script performs GTA L and K Pairwise Compares for user specified pairs of Experiments
|
|
#
|
|
# Updated 240724 Bryan C Roessler to improve file operations and portability
|
|
# NOTE: The two required arguments are the same and now there are two optional arguments
|
|
# 1. Exp1
|
|
# 2. Exp2
|
|
# 3. StudyInfo.csv file
|
|
# 4. Output Directory
|
|
|
|
library("ggplot2")
|
|
library("plotly")
|
|
library("htmlwidgets")
|
|
library("extrafont")
|
|
library("grid")
|
|
library("ggthemes")
|
|
|
|
args <- commandArgs(TRUE)
|
|
exp_name <- args[1]
|
|
exp_name2 <- args[2]
|
|
|
|
if (length(args) >= 3) {
|
|
study_info_file <- args[3]
|
|
} else {
|
|
study_info_file <- "StudyInfo.csv"
|
|
}
|
|
|
|
if (length(args) > 4) {
|
|
output_dir <- args[4]
|
|
} else {
|
|
output_dir <- "gta"
|
|
}
|
|
|
|
# if (length(args) > 3) {
|
|
# sgd_terms_file <- args[3]
|
|
# } else {
|
|
# sgd_terms_file <- "../Code/go_terms.tab"
|
|
# }
|
|
|
|
# if (length(args) > 4) {
|
|
# sgd_features_file <- args[4]
|
|
# } else {
|
|
# sgd_features_file <- "../Code/gene_association.sgd"
|
|
# }
|
|
|
|
expNumber1<- as.numeric(sub("^.*?(\\d+)$", "\\1", exp_name))
|
|
expNumber2<- as.numeric(sub("^.*?(\\d+)$", "\\1", exp_name2))
|
|
Labels<- read.csv(file=study_info_file,stringsAsFactors = FALSE)
|
|
Name1 <- Labels[expNumber1,2]
|
|
Name2 <- Labels[expNumber2,2]
|
|
|
|
go_terms_file <- "Average_GOTerms_All.csv"
|
|
input_file1 <- file.path(output_dir,exp_name,go_terms_file)
|
|
input_file2 <- file.path(output_dir,exp_name2,go_terms_file)
|
|
|
|
pairDirL= file.path(output_dir,"PairwiseCompareL_",exp_name,"-",exp_name2)
|
|
pairDirK= file.path(output_dir,"PairwiseCompareK_",exp_name,"-",exp_name2)
|
|
outPathGTAcompare= file.path(output_dir,"PairwiseCompareL")
|
|
outPathGTAcompare[2]= file.path(output_dir,"PairwiseCompareK")
|
|
#dir.create(outPathGTAcompare[1])
|
|
dir.create(pairDirL) #(outPathGTAcompare[1])
|
|
dir.create(pairDirK) #(outPathGTAcompare[2])
|
|
|
|
###########BEGIN PAIRWISE L-----LLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLL
|
|
outputpath <- pairDirL #outPathGTAcompare[1] #Args[5]
|
|
#outputPlotly <- "../GTAresults/PairwiseCompareL/" #"/GTAresults/PairwiseCompareL/"
|
|
print("39")
|
|
#theme elements for plots
|
|
theme_Publication <- function(base_size=14, base_family="sans") {
|
|
(theme_foundation(base_size=base_size, base_family=base_family)
|
|
+ theme(plot.title = element_text(face = "bold",
|
|
size = rel(1.2), hjust = 0.5),
|
|
text = element_text(),
|
|
panel.background = element_rect(colour = NA),
|
|
plot.background = element_rect(colour = NA),
|
|
panel.border = element_rect(colour = NA),
|
|
axis.title = element_text(face = "bold",size = rel(1)),
|
|
axis.title.y = element_text(angle=90,vjust =2),
|
|
axis.title.x = element_text(vjust = -0.2),
|
|
axis.text = element_text(),
|
|
axis.line = element_line(colour="black"),
|
|
axis.ticks = element_line(),
|
|
panel.grid.major = element_line(colour="#f0f0f0"),
|
|
panel.grid.minor = element_blank(),
|
|
legend.key = element_rect(colour = NA),
|
|
legend.position = "bottom",
|
|
legend.direction = "horizontal",
|
|
legend.key.size= unit(0.2, "cm"),
|
|
legend.spacing = unit(0, "cm"),
|
|
legend.title = element_text(face="italic"),
|
|
plot.margin=unit(c(10,5,5,5),"mm"),
|
|
strip.background=element_rect(colour="#f0f0f0",fill="#f0f0f0"),
|
|
strip.text = element_text(face="bold")
|
|
))
|
|
|
|
}
|
|
|
|
scale_fill_Publication <- function(...){
|
|
library(scales)
|
|
discrete_scale("fill","Publication",manual_pal(values = c("#386cb0","#fdb462","#7fc97f","#ef3b2c","#662506","#a6cee3","#fb9a99","#984ea3","#ffff33")), ...)
|
|
}
|
|
|
|
scale_colour_Publication <- function(...){
|
|
discrete_scale("colour","Publication",manual_pal(values = c("#386cb0","#fdb462","#7fc97f","#ef3b2c","#662506","#a6cee3","#fb9a99","#984ea3","#ffff33")), ...)
|
|
}
|
|
|
|
theme_Publication_legend_right <- function(base_size=14, base_family="sans") {
|
|
(theme_foundation(base_size=base_size, base_family=base_family)
|
|
+ theme(plot.title = element_text(face = "bold",
|
|
size = rel(1.2), hjust = 0.5),
|
|
text = element_text(),
|
|
panel.background = element_rect(colour = NA),
|
|
plot.background = element_rect(colour = NA),
|
|
panel.border = element_rect(colour = NA),
|
|
axis.title = element_text(face = "bold",size = rel(1)),
|
|
axis.title.y = element_text(angle=90,vjust =2),
|
|
axis.title.x = element_text(vjust = -0.2),
|
|
axis.text = element_text(),
|
|
axis.line = element_line(colour="black"),
|
|
axis.ticks = element_line(),
|
|
panel.grid.major = element_line(colour="#f0f0f0"),
|
|
panel.grid.minor = element_blank(),
|
|
legend.key = element_rect(colour = NA),
|
|
legend.position = "right",
|
|
legend.direction = "vertical",
|
|
legend.key.size= unit(0.5, "cm"),
|
|
legend.spacing = unit(0, "cm"),
|
|
legend.title = element_text(face="italic"),
|
|
plot.margin=unit(c(10,5,5,5),"mm"),
|
|
strip.background=element_rect(colour="#f0f0f0",fill="#f0f0f0"),
|
|
strip.text = element_text(face="bold")
|
|
))
|
|
}
|
|
|
|
scale_fill_Publication <- function(...){
|
|
discrete_scale("fill","Publication",manual_pal(values = c("#386cb0","#fdb462","#7fc97f","#ef3b2c","#662506","#a6cee3","#fb9a99","#984ea3","#ffff33")), ...)
|
|
}
|
|
|
|
scale_colour_Publication <- function(...){
|
|
discrete_scale("colour","Publication",manual_pal(values = c("#386cb0","#fdb462","#7fc97f","#ef3b2c","#662506","#a6cee3","#fb9a99","#984ea3","#ffff33")), ...)
|
|
}
|
|
|
|
X1 <- read.csv(file = input_file1,stringsAsFactors=FALSE,header = TRUE)
|
|
X2 <- read.csv(file = input_file2,stringsAsFactors=FALSE,header = TRUE)
|
|
print(117)
|
|
|
|
X <- merge(X1,X2,by ="Term_Avg",all=TRUE,suffixes = c("_X1","_X2"))
|
|
gg <- ggplot(data = X,aes(x=Z_lm_L_Avg_X1,y=Z_lm_L_Avg_X2,color=Ontology_Avg_X1,Term=Term_Avg,Genes=Genes_Avg_X1,NumGenes=NumGenes_Avg_X1,AllPossibleGenes=AllPossibleGenes_Avg_X1,SD_1=Z_lm_L_SD_X1,SD_2=Z_lm_L_SD_X2)) +
|
|
xlab(paste("GO Term Avg lm Z for ",Name1,sep="")) +
|
|
geom_rect(aes(xmin=-2,xmax=2,ymin=-2,ymax=2),color="grey20",size=0.25,alpha=0.1,inherit.aes = FALSE,fill=NA) + geom_point(shape=3) +
|
|
ylab(paste("GO Term Avg lm Z for ",Name2,sep="")) + ggtitle(paste("Comparing Average GO Term Z lm for ",Name1," vs. ",Name2,sep = "")) +
|
|
theme_Publication_legend_right()
|
|
pdf(paste(outputpath,"/","Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_ByOntology.pdf",sep=""),width = 12,height = 8)
|
|
gg
|
|
dev.off()
|
|
pgg <- ggplotly(gg)
|
|
#pgg
|
|
print("130")
|
|
fname <- paste("Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_byOntology.html",sep="")
|
|
print(fname)
|
|
htmlwidgets::saveWidget(pgg, file.path(getwd(),fname))
|
|
file.rename(from = file.path(getwd(),fname), to = file.path(pairDirL,fname))
|
|
|
|
|
|
#ID aggravators and alleviators, regardless of whether they meet 2SD threshold
|
|
X1_Specific_Aggravators <- X[which(X$Z_lm_L_Avg_X1 >= 2 & X$Z_lm_L_Avg_X2 < 2),]
|
|
X1_Specific_Alleviators <- X[which(X$Z_lm_L_Avg_X1 <= -2 & X$Z_lm_L_Avg_X2 > -2),]
|
|
|
|
X2_Specific_Aggravators <- X[which(X$Z_lm_L_Avg_X2 >= 2 & X$Z_lm_L_Avg_X1 < 2),]
|
|
X2_Specific_Alleviators <- X[which(X$Z_lm_L_Avg_X2 <= -2 & X$Z_lm_L_Avg_X1 > -2),]
|
|
|
|
Overlap_Aggravators <- X[which(X$Z_lm_L_Avg_X1 >= 2 & X$Z_lm_L_Avg_X2 >= 2),]
|
|
Overlap_Alleviators <- X[which(X$Z_lm_L_Avg_X1 <= -2 & X$Z_lm_L_Avg_X2 <= -2),]
|
|
|
|
X2_Specific_Aggravators_X1_Alleviatiors <- X[which(X$Z_lm_L_Avg_X2 >= 2 & X$Z_lm_L_Avg_X1 <= -2),]
|
|
X2_Specific_Alleviators_X1_Aggravators <- X[which(X$Z_lm_L_Avg_X2 <= -2 & X$Z_lm_L_Avg_X1 >= 2),]
|
|
print("155")
|
|
X$Overlap_Avg <- NA
|
|
|
|
try(X[X$Term_Avg %in% X1_Specific_Aggravators$Term_Avg,]$Overlap_Avg <- paste(Name1,"Specific_Deletion_Enhancers",sep="_"))
|
|
try(X[X$Term_Avg %in% X1_Specific_Alleviators$Term_Avg,]$Overlap_Avg <- paste(Name1,"Specific_Deletion_Suppresors",sep="_"))
|
|
try(X[X$Term_Avg %in% X2_Specific_Aggravators$Term_Avg,]$Overlap_Avg <- paste(Name2,"Specific_Deletion_Enhancers",sep="_"))
|
|
try(X[X$Term_Avg %in% X2_Specific_Alleviators$Term_Avg,]$Overlap_Avg <- paste(Name2,"Specific_Deletion_Suppressors",sep="_"))
|
|
try(X[X$Term_Avg %in% Overlap_Aggravators$Term_Avg,]$Overlap_Avg <- "Overlapping_Deletion_Enhancers")
|
|
try(X[X$Term_Avg %in% Overlap_Alleviators$Term_Avg,]$Overlap_Avg <- "Overlapping_Deletion_Suppressors")
|
|
try(X[X$Term_Avg %in% X2_Specific_Aggravators_X1_Alleviatiors$Term_Avg,]$Overlap_Avg <- paste(Name2,"Deletion_Enhancers",Name1,"Deletion_Suppressors",sep="_"))
|
|
try(X[X$Term_Avg %in% X2_Specific_Alleviators_X1_Aggravators$Term_Avg,]$Overlap_Avg <- paste(Name2,"Deletion_Suppressors",Name1,"Deletion_Enhancers",sep="_"))
|
|
|
|
|
|
gg <- ggplot(data = X,aes(x=Z_lm_L_Avg_X1,y=Z_lm_L_Avg_X2,color=Overlap_Avg,Term=Term_Avg,Genes=Genes_Avg_X1,NumGenes=NumGenes_Avg_X1,AllPossibleGenes=AllPossibleGenes_Avg_X1,SD_1=Z_lm_L_SD_X1,SD_2=Z_lm_L_SD_X2)) +
|
|
xlab(paste("GO Term Avg lm Z for ",Name1,sep="")) +
|
|
geom_rect(aes(xmin=-2,xmax=2,ymin=-2,ymax=2),color="grey20",size=0.25,alpha=0.1,inherit.aes = FALSE,fill=NA) + geom_point(shape=3) +
|
|
ylab(paste("GO Term Avg lm Z for ",Name2,sep="")) + ggtitle(paste("Comparing Average GO Term Z lm for ",Name1," vs. ",Name2,sep="")) +
|
|
theme_Publication_legend_right()
|
|
pdf(paste(outputpath,"/","Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_ByOverlap.pdf",sep=""),width = 12,height = 8)
|
|
gg
|
|
dev.off()
|
|
pgg <- ggplotly(gg)
|
|
|
|
#pgg
|
|
print("#2-174")
|
|
#2
|
|
fname <- paste("/Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_byOverlap.html",sep="")
|
|
print(fname)
|
|
htmlwidgets::saveWidget(pgg, file.path(getwd(),fname))
|
|
file.rename(from = file.path(getwd(),fname), to = file.path(pairDirL,fname))
|
|
|
|
x_rem2_gene <- X[X$NumGenes_Avg_X1 >= 2 & X$NumGenes_Avg_X2 >= 2,]
|
|
#3
|
|
gg <- ggplot(data = x_rem2_gene,aes(x=Z_lm_L_Avg_X1,y=Z_lm_L_Avg_X2,color=Overlap_Avg,Term=Term_Avg,Genes=Genes_Avg_X1,NumGenes=NumGenes_Avg_X1,AllPossibleGenes=AllPossibleGenes_Avg_X1,SD_1=Z_lm_L_SD_X1,SD_2=Z_lm_L_SD_X2)) +
|
|
xlab(paste("GO Term Avg lm Z for ",Name1,sep="")) +
|
|
geom_rect(aes(xmin=-2,xmax=2,ymin=-2,ymax=2),color="grey20",size=0.25,alpha=0.1,inherit.aes = FALSE,fill=NA) + geom_point(shape=3) +
|
|
ylab(paste("GO Term Avg lm Z for ",Name2,sep="")) + ggtitle(paste("Comparing Average GO Term Z lm for ",Name1," vs. ",Name2,sep="")) +
|
|
theme_Publication_legend_right()
|
|
pdf(paste(outputpath,"/","Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_ByOverlap_above2genes.pdf",sep=""),width = 12,height = 8)
|
|
gg
|
|
dev.off()
|
|
pgg <- ggplotly(gg)
|
|
#pgg
|
|
print("#3")
|
|
fname <- paste("Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_byOverlap_above2genes.html",sep="")
|
|
print(fname)
|
|
htmlwidgets::saveWidget(pgg, file.path(getwd(),fname))
|
|
file.rename(from = file.path(getwd(),fname), to = file.path(pairDirL,fname))
|
|
|
|
#4
|
|
X_overlap_nothresold <- X[!(is.na(X$Overlap_Avg)),]
|
|
gg <- ggplot(data = X_overlap_nothresold,aes(x=Z_lm_L_Avg_X1,y=Z_lm_L_Avg_X2,color=Overlap_Avg,Term=Term_Avg,Genes=Genes_Avg_X1,NumGenes=NumGenes_Avg_X1,AllPossibleGenes=AllPossibleGenes_Avg_X1,SD_1=Z_lm_L_SD_X1,SD_2=Z_lm_L_SD_X2)) +
|
|
xlab(paste("GO Term Avg lm Z for ",Name1,sep="")) +
|
|
geom_rect(aes(xmin=-2,xmax=2,ymin=-2,ymax=2),color="grey20",size=0.25,alpha=0.1,inherit.aes = FALSE,fill=NA) + geom_point(shape=3) +
|
|
ylab(paste("GO Term Avg lm Z for ",Name2,sep="")) + ggtitle(paste("Comparing Average GO Term Z lm for ",Name1," vs. ",Name2,sep = "")) +
|
|
theme_Publication_legend_right()
|
|
pdf(paste(outputpath,"/","Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_Above2SD_ByOverlap.pdf",sep=""),width = 12,height = 8)
|
|
gg
|
|
dev.off()
|
|
pgg <- ggplotly(gg)
|
|
#pgg
|
|
print("#4")
|
|
fname <- paste("Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_Above2SD_ByOverlap.html",sep="")
|
|
print(fname)
|
|
htmlwidgets::saveWidget(pgg, file.path(getwd(),fname))
|
|
file.rename(from = file.path(getwd(),fname), to = file.path(pairDirL,fname))
|
|
|
|
|
|
#only output GTA terms where average score is still above 2 after subtracting the SD
|
|
#Z1 will ID aggravators, Z2 alleviators
|
|
Z1 <- X
|
|
Z1$L_Subtract_SD_X1 <- Z1$Z_lm_L_Avg_X1 - Z1$Z_lm_L_SD_X1
|
|
Z1$L_Subtract_SD_X2 <- Z1$Z_lm_L_Avg_X2 - Z1$Z_lm_L_SD_X2
|
|
|
|
Z2 <- X
|
|
Z2$L_Subtract_SD_X1 <- Z1$Z_lm_L_Avg_X1 + Z1$Z_lm_L_SD_X1
|
|
Z2$L_Subtract_SD_X2 <- Z1$Z_lm_L_Avg_X2 + Z1$Z_lm_L_SD_X2
|
|
|
|
X1_Specific_Aggravators2 <- Z1[which(Z1$L_Subtract_SD_X1 >= 2 & Z1$L_Subtract_SD_X2 < 2),]
|
|
X1_Specific_Alleviators2 <- Z2[which(Z2$L_Subtract_SD_X1 <= -2 & Z2$L_Subtract_SD_X2 > -2),]
|
|
|
|
X2_Specific_Aggravators2 <- Z1[which(Z1$L_Subtract_SD_X2 >= 2 & Z1$L_Subtract_SD_X1 < 2),]
|
|
X2_Specific_Alleviators2 <- Z2[which(Z2$L_Subtract_SD_X2 <= -2 & Z2$L_Subtract_SD_X1 > -2),]
|
|
|
|
Overlap_Aggravators2 <- Z1[which(Z1$L_Subtract_SD_X1 >= 2 & Z1$L_Subtract_SD_X2 >= 2),]
|
|
Overlap_Alleviators2 <- Z2[which(Z2$L_Subtract_SD_X2 <= -2 & Z2$L_Subtract_SD_X1 <= -2),]
|
|
|
|
X2_Specific_Aggravators2_X1_Alleviatiors2 <- Z1[which(Z1$L_Subtract_SD_X2 >= 2 & Z2$L_Subtract_SD_X1 <= -2),]
|
|
X2_Specific_Alleviators2_X1_Aggravators2 <- Z2[which(Z2$L_Subtract_SD_X2 <= -2 & Z1$L_Subtract_SD_X1 >= 2),]
|
|
|
|
X$Overlap <- NA
|
|
|
|
try(X[X$Term_Avg %in% X1_Specific_Aggravators2$Term_Avg,]$Overlap <- paste(Name1,"Specific_Deletion_Enhancers",sep="_"))
|
|
try(X[X$Term_Avg %in% X1_Specific_Alleviators2$Term_Avg,]$Overlap <- paste(Name1,"Specific_Deletion_Suppresors",sep="_"))
|
|
try(X[X$Term_Avg %in% X2_Specific_Aggravators2$Term_Avg,]$Overlap <- paste(Name2,"Specific_Deletion_Enhancers",sep="_"))
|
|
try(X[X$Term_Avg %in% X2_Specific_Alleviators2$Term_Avg,]$Overlap <- paste(Name2,"Specific_Deletion_Suppressors",sep="_"))
|
|
try(X[X$Term_Avg %in% Overlap_Aggravators2$Term_Avg,]$Overlap <- "Overlapping_Deletion_Enhancers")
|
|
try(X[X$Term_Avg %in% Overlap_Alleviators2$Term_Avg,]$Overlap <- "Overlapping_Deletion_Suppressors")
|
|
try(X[X$Term_Avg %in% X2_Specific_Aggravators2_X1_Alleviatiors2$Term_Avg,]$Overlap <- paste(Name2,"Deletion_Enhancers",Name1,"Deletion_Suppressors",sep="_"))
|
|
try(X[X$Term_Avg %in% X2_Specific_Alleviators2_X1_Aggravators2$Term_Avg,]$Overlap <- paste(Name2,"Deletion_Suppressors",Name1,"Deletion_Enhancers",sep="_"))
|
|
|
|
#5
|
|
X_abovethreshold <- X[!(is.na(X$Overlap)),]
|
|
|
|
gg <- ggplot(data = X_abovethreshold,aes(x=Z_lm_L_Avg_X1,y=Z_lm_L_Avg_X2,color=Overlap,Term=Term_Avg,Genes=Genes_Avg_X1,NumGenes=NumGenes_Avg_X1,AllPossibleGenes=AllPossibleGenes_Avg_X1,SD_1=Z_lm_L_SD_X1,SD_2=Z_lm_L_SD_X2)) +
|
|
xlab(paste("GO Term Avg lm Z for ",Name1,sep="")) +
|
|
geom_rect(aes(xmin=-2,xmax=2,ymin=-2,ymax=2),color="grey20",size=0.25,alpha=0.1,inherit.aes = FALSE,fill=NA) + geom_point(shape=3) +
|
|
ylab(paste("GO Term Avg lm Z for ",Name2,sep="")) + ggtitle(paste("Comparing Average GO Term Z lm for ",Name1," vs. ",Name2,sep="")) +
|
|
theme_Publication_legend_right()
|
|
pdf(paste(outputpath,"/","Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_ByOverlap_AboveThreshold.pdf",sep=""),width = 12,height = 8)
|
|
gg
|
|
dev.off()
|
|
pgg <- ggplotly(gg)
|
|
#pgg
|
|
print("#5")
|
|
fname <- paste("Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_ByOverlap_AboveThreshold.html",sep="")
|
|
print(fname)
|
|
htmlwidgets::saveWidget(pgg, file.path(getwd(),fname))
|
|
file.rename(from = file.path(getwd(),fname), to = file.path(pairDirL,fname))
|
|
#6
|
|
gg <- ggplot(data = X_abovethreshold,aes(x=Z_lm_L_Avg_X1,y=Z_lm_L_Avg_X2,color=Overlap,Term=Term_Avg,Genes=Genes_Avg_X1,NumGenes=NumGenes_Avg_X1,AllPossibleGenes=AllPossibleGenes_Avg_X1,SD_1=Z_lm_L_SD_X1,SD_2=Z_lm_L_SD_X2)) +
|
|
xlab(paste("GO Term Avg lm Z for ",Name1,sep="")) +
|
|
geom_text(aes(label=Term_Avg),nudge_y = 0.25,size=2) +
|
|
geom_rect(aes(xmin=-2,xmax=2,ymin=-2,ymax=2),color="grey20",size=0.25,alpha=0.1,inherit.aes = FALSE,fill=NA) + geom_point(shape=3,size=3) +
|
|
ylab(paste("GO Term Avg lm Z for ",Name2,sep="")) + ggtitle(paste("Comparing Average GO Term Z lm for ",Name1," vs. ",Name2,sep="")) +
|
|
theme_Publication_legend_right()
|
|
pdf(paste(outputpath,"/","Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_ByOverlap_AboveThreshold_names.pdf",sep=""),width = 20,height = 20)
|
|
gg
|
|
dev.off()
|
|
pgg <- ggplotly(gg)
|
|
#pgg
|
|
print("#6")
|
|
fname <- paste("Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_ByOverlap_AboveThreshold_names.html",sep="")
|
|
print(fname)
|
|
htmlwidgets::saveWidget(pgg, file.path(getwd(),fname))
|
|
file.rename(from = file.path(getwd(),fname), to = file.path(pairDirL,fname))
|
|
|
|
X_abovethreshold$X1_Rank <- NA
|
|
X_abovethreshold$X1_Rank <- rank(-X_abovethreshold$Z_lm_L_Avg_X1,ties.method = "random")
|
|
X_abovethreshold$X2_Rank <- NA
|
|
X_abovethreshold$X2_Rank <- rank(-X_abovethreshold$Z_lm_L_Avg_X2,ties.method = "random")
|
|
|
|
#7
|
|
gg <- ggplot(data = X_abovethreshold,aes(x=Z_lm_L_Avg_X1,y=Z_lm_L_Avg_X2,color=Overlap,Term=Term_Avg,Genes=Genes_Avg_X1,NumGenes=NumGenes_Avg_X1,AllPossibleGenes=AllPossibleGenes_Avg_X1,SD_1=Z_lm_L_SD_X1,SD_2=Z_lm_L_SD_X2)) +
|
|
xlab(paste("GO Term Avg lm Z for ",Name1,sep="")) +
|
|
geom_text(aes(label=X1_Rank),nudge_y = 0.25,size=4) +
|
|
geom_rect(aes(xmin=-2,xmax=2,ymin=-2,ymax=2),color="grey20",size=0.25,alpha=0.1,inherit.aes = FALSE,fill=NA) + geom_point(shape=3,size=3) +
|
|
ylab(paste("GO Term Avg lm Z for ",Name2,sep="")) + ggtitle(paste("Comparing Average GO Term Z lm for ",Name1," vs. ",Name2,sep="")) +
|
|
theme_Publication_legend_right()
|
|
pdf(paste(outputpath,"/","Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_ByOverlap_AboveThreshold_numberedX1.pdf",sep=""),width = 15,height = 15)
|
|
gg
|
|
dev.off()
|
|
pgg <- ggplotly(gg)
|
|
#pgg
|
|
print("#7")
|
|
fname <- paste("Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_ByOverlap_AboveThreshold_numberedX1.html",sep="")
|
|
print(fname)
|
|
htmlwidgets::saveWidget(pgg, file.path(getwd(),fname))
|
|
file.rename(from = file.path(getwd(),fname), to = file.path(pairDirL,fname))
|
|
|
|
#8
|
|
gg <- ggplot(data = X_abovethreshold,aes(x=Z_lm_L_Avg_X1,y=Z_lm_L_Avg_X2,color=Overlap,Term=Term_Avg,Genes=Genes_Avg_X1,NumGenes=NumGenes_Avg_X1,AllPossibleGenes=AllPossibleGenes_Avg_X1,SD_1=Z_lm_L_SD_X1,SD_2=Z_lm_L_SD_X2)) +
|
|
xlab(paste("GO Term Avg lm Z for ",Name1,sep="")) +
|
|
geom_text(aes(label=X2_Rank),nudge_y = 0.25,size=4) +
|
|
geom_rect(aes(xmin=-2,xmax=2,ymin=-2,ymax=2),color="grey20",size=0.25,alpha=0.1,inherit.aes = FALSE,fill=NA) + geom_point(shape=3,size=3) +
|
|
ylab(paste("GO Term Avg lm Z for ",Name2,sep="")) + ggtitle(paste("Comparing Average GO Term Z lm for ",Name1," vs. ",Name2,sep="")) +
|
|
theme_Publication_legend_right()
|
|
pdf(paste(outputpath,"/","Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_ByOverlap_AboveThreshold_numberedX2.pdf",sep=""),width = 15,height = 15)
|
|
gg
|
|
dev.off()
|
|
pgg <- ggplotly(gg)
|
|
#pgg
|
|
print("#8")
|
|
fname <- paste("Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_ByOverlap_AboveThreshold_numberedX2.html",sep="")
|
|
print(fname)
|
|
htmlwidgets::saveWidget(pgg, file.path(getwd(),fname))
|
|
file.rename(from = file.path(getwd(),fname), to = file.path(pairDirL,fname))
|
|
|
|
print("write csv files L")
|
|
write.csv(x=X,file = paste(outputpath,"/All_GTA_Avg_Scores_",Name1,"_vs_",Name2,".csv",sep=""),row.names = FALSE)
|
|
write.csv(x=X_abovethreshold,file = paste(getwd(),"/",outputpath,"/","AboveThreshold_GTA_Avg_Scores_",Name1,"_vs_",Name2,".csv",sep=""),row.names = FALSE)
|
|
#End of L GTA Pairwise Compare
|
|
|
|
|
|
|
|
###########BEGIN PAIRWISE K-----KKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKK
|
|
#define the output path (as fourth argument from Rscript)
|
|
outputpath <- pairDirK #outPathGTAcompare[2] #Args[5]
|
|
|
|
#theme elements for plots
|
|
theme_Publication <- function(base_size=14, base_family="sans") {
|
|
(theme_foundation(base_size=base_size, base_family=base_family)
|
|
+ theme(plot.title = element_text(face = "bold",
|
|
size = rel(1.2), hjust = 0.5),
|
|
text = element_text(),
|
|
panel.background = element_rect(colour = NA),
|
|
plot.background = element_rect(colour = NA),
|
|
panel.border = element_rect(colour = NA),
|
|
axis.title = element_text(face = "bold",size = rel(1)),
|
|
axis.title.y = element_text(angle=90,vjust =2),
|
|
axis.title.x = element_text(vjust = -0.2),
|
|
axis.text = element_text(),
|
|
axis.line = element_line(colour="black"),
|
|
axis.ticks = element_line(),
|
|
panel.grid.major = element_line(colour="#f0f0f0"),
|
|
panel.grid.minor = element_blank(),
|
|
legend.key = element_rect(colour = NA),
|
|
legend.position = "bottom",
|
|
legend.direction = "horizontal",
|
|
legend.key.size= unit(0.2, "cm"),
|
|
legend.spacing = unit(0, "cm"),
|
|
legend.title = element_text(face="italic"),
|
|
plot.margin=unit(c(10,5,5,5),"mm"),
|
|
strip.background=element_rect(colour="#f0f0f0",fill="#f0f0f0"),
|
|
strip.text = element_text(face="bold")
|
|
))
|
|
}
|
|
|
|
scale_fill_Publication <- function(...){
|
|
library(scales)
|
|
discrete_scale("fill","Publication",manual_pal(values = c("#386cb0","#fdb462","#7fc97f","#ef3b2c","#662506","#a6cee3","#fb9a99","#984ea3","#ffff33")), ...)
|
|
}
|
|
scale_colour_Publication <- function(...){
|
|
discrete_scale("colour","Publication",manual_pal(values = c("#386cb0","#fdb462","#7fc97f","#ef3b2c","#662506","#a6cee3","#fb9a99","#984ea3","#ffff33")), ...)
|
|
}
|
|
|
|
theme_Publication_legend_right <- function(base_size=14, base_family="sans") {
|
|
(theme_foundation(base_size=base_size, base_family=base_family)
|
|
+ theme(plot.title = element_text(face = "bold",
|
|
size = rel(1.2), hjust = 0.5),
|
|
text = element_text(),
|
|
panel.background = element_rect(colour = NA),
|
|
plot.background = element_rect(colour = NA),
|
|
panel.border = element_rect(colour = NA),
|
|
axis.title = element_text(face = "bold",size = rel(1)),
|
|
axis.title.y = element_text(angle=90,vjust =2),
|
|
axis.title.x = element_text(vjust = -0.2),
|
|
axis.text = element_text(),
|
|
axis.line = element_line(colour="black"),
|
|
axis.ticks = element_line(),
|
|
panel.grid.major = element_line(colour="#f0f0f0"),
|
|
panel.grid.minor = element_blank(),
|
|
legend.key = element_rect(colour = NA),
|
|
legend.position = "right",
|
|
legend.direction = "vertical",
|
|
legend.key.size= unit(0.5, "cm"),
|
|
legend.spacing = unit(0, "cm"),
|
|
legend.title = element_text(face="italic"),
|
|
plot.margin=unit(c(10,5,5,5),"mm"),
|
|
strip.background=element_rect(colour="#f0f0f0",fill="#f0f0f0"),
|
|
strip.text = element_text(face="bold")
|
|
))
|
|
}
|
|
|
|
scale_fill_Publication <- function(...){
|
|
discrete_scale("fill","Publication",manual_pal(values = c("#386cb0","#fdb462","#7fc97f","#ef3b2c","#662506","#a6cee3","#fb9a99","#984ea3","#ffff33")), ...)
|
|
}
|
|
scale_colour_Publication <- function(...){
|
|
discrete_scale("colour","Publication",manual_pal(values = c("#386cb0","#fdb462","#7fc97f","#ef3b2c","#662506","#a6cee3","#fb9a99","#984ea3","#ffff33")), ...)
|
|
}
|
|
|
|
X1 <- read.csv(file = input_file1,stringsAsFactors=FALSE,header = TRUE)
|
|
|
|
X2 <- read.csv(file = input_file2,stringsAsFactors=FALSE,header = TRUE)
|
|
|
|
#1
|
|
X <- merge(X1,X2,by ="Term_Avg",all=TRUE,suffixes = c("_X1","_X2"))
|
|
|
|
gg <- ggplot(data = X,aes(x=Z_lm_K_Avg_X1,y=Z_lm_K_Avg_X2,color=Ontology_Avg_X1,Term=Term_Avg,Genes=Genes_Avg_X1,NumGenes=NumGenes_Avg_X1,AllPossibleGenes=AllPossibleGenes_Avg_X1,SD_1=Z_lm_K_SD_X1,SD_2=Z_lm_K_SD_X2)) +
|
|
xlab(paste("GO Term Avg lm Z for ",Name1,sep="")) +
|
|
geom_rect(aes(xmin=-2,xmax=2,ymin=-2,ymax=2),color="grey20",size=0.25,alpha=0.1,inherit.aes = FALSE,fill=NA) + geom_point(shape=3) +
|
|
ylab(paste("GO Term Avg lm Z for ",Name2,sep="")) + ggtitle(paste("Comparing Average GO Term Z lm for ",Name1," vs. ",Name2,sep = "")) +
|
|
theme_Publication_legend_right()
|
|
pdf(paste(getwd(),"/",outputpath,"/","Scatter_lm_GTF_Averages_",Name1,"_vs_",Name2,"_All_ByOntology.pdf",sep=""),width = 12,height = 8)
|
|
gg
|
|
dev.off()
|
|
pgg <- ggplotly(gg)
|
|
#pgg
|
|
print("127")
|
|
fname <- paste("Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_byOntology.html",sep="")
|
|
print(fname)
|
|
htmlwidgets::saveWidget(pgg, file.path(getwd(),fname))
|
|
file.rename(from = file.path(getwd(),fname), to = file.path(pairDirK,fname))
|
|
|
|
#2
|
|
#ID aggravators and alleviators, regardless of whether they meet 2SD threshold
|
|
X1_Specific_Aggravators <- X[which(X$Z_lm_K_Avg_X1 >= 2 & X$Z_lm_K_Avg_X2 < 2),]
|
|
X1_Specific_Alleviators <- X[which(X$Z_lm_K_Avg_X1 <= -2 & X$Z_lm_K_Avg_X2 > -2),]
|
|
|
|
X2_Specific_Aggravators <- X[which(X$Z_lm_K_Avg_X2 >= 2 & X$Z_lm_K_Avg_X1 < 2),]
|
|
X2_Specific_Alleviators <- X[which(X$Z_lm_K_Avg_X2 <= -2 & X$Z_lm_K_Avg_X1 > -2),]
|
|
|
|
Overlap_Aggravators <- X[which(X$Z_lm_K_Avg_X1 >= 2 & X$Z_lm_K_Avg_X2 >= 2),]
|
|
Overlap_Alleviators <- X[which(X$Z_lm_K_Avg_X1 <= -2 & X$Z_lm_K_Avg_X2 <= -2),]
|
|
|
|
X2_Specific_Aggravators_X1_Alleviatiors <- X[which(X$Z_lm_K_Avg_X2 >= 2 & X$Z_lm_K_Avg_X1 <= -2),]
|
|
X2_Specific_Alleviators_X1_Aggravators <- X[which(X$Z_lm_K_Avg_X2 <= -2 & X$Z_lm_K_Avg_X1 >= 2),]
|
|
|
|
X$Overlap_Avg <- NA
|
|
|
|
try(X[X$Term_Avg %in% X1_Specific_Aggravators$Term_Avg,]$Overlap_Avg <- paste(Name1,"Specific_Deletion_Suppressors",sep="_"))
|
|
try(X[X$Term_Avg %in% X1_Specific_Alleviators$Term_Avg,]$Overlap_Avg <- paste(Name1,"Specific_Deletion_Enhancers",sep="_"))
|
|
try(X[X$Term_Avg %in% X2_Specific_Aggravators$Term_Avg,]$Overlap_Avg <- paste(Name2,"Specific_Deletion_Suppressors",sep="_"))
|
|
try(X[X$Term_Avg %in% X2_Specific_Alleviators$Term_Avg,]$Overlap_Avg <- paste(Name2,"Specific_Deletion_Enhancers",sep="_"))
|
|
try(X[X$Term_Avg %in% Overlap_Aggravators$Term_Avg,]$Overlap_Avg <- "Overlapping_Deletion_Suppressors")
|
|
try(X[X$Term_Avg %in% Overlap_Alleviators$Term_Avg,]$Overlap_Avg <- "Overlapping_Deletion_Enhancers")
|
|
try(X[X$Term_Avg %in% X2_Specific_Aggravators_X1_Alleviatiors$Term_Avg,]$Overlap_Avg <- paste(Name2,"Deletion_Suppressors",Name1,"Deletion_Enhancers",sep="_"))
|
|
try(X[X$Term_Avg %in% X2_Specific_Alleviators_X1_Aggravators$Term_Avg,]$Overlap_Avg <- paste(Name2,"Deletion_Enhancers",Name1,"Deletion_Suppressors",sep="_"))
|
|
|
|
|
|
plotly_path <- paste(getwd(),"/",outputpath,"/","Scatter_lm_GTF_Averages_",Name1,"_vs_",Name2,"_All_byOverlap.html",sep="")
|
|
gg <- ggplot(data = X,aes(x=Z_lm_K_Avg_X1,y=Z_lm_K_Avg_X2,color=Overlap_Avg,Term=Term_Avg,Genes=Genes_Avg_X1,NumGenes=NumGenes_Avg_X1,AllPossibleGenes=AllPossibleGenes_Avg_X1,SD_1=Z_lm_K_SD_X1,SD_2=Z_lm_K_SD_X2)) +
|
|
xlab(paste("GO Term Avg lm Z for ",Name1,sep="")) +
|
|
geom_rect(aes(xmin=-2,xmax=2,ymin=-2,ymax=2),color="grey20",size=0.25,alpha=0.1,inherit.aes = FALSE,fill=NA) + geom_point(shape=3) +
|
|
ylab(paste("GO Term Avg lm Z for ",Name2,sep="")) + ggtitle(paste("Comparing Average GO Term Z lm for ",Name1," vs. ",Name2,sep="")) +
|
|
theme_Publication_legend_right()
|
|
pdf(paste(getwd(),"/",outputpath,"/","Scatter_lm_GTF_Averages_",Name1,"_vs_",Name2,"_All_ByOverlap.pdf",sep=""),width = 12,height = 8)
|
|
gg
|
|
dev.off()
|
|
pgg <- ggplotly(gg)
|
|
#pgg
|
|
print("#2-170")
|
|
#2
|
|
fname <- paste("/Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_byOverlap.html",sep="")
|
|
print(fname)
|
|
htmlwidgets::saveWidget(pgg, file.path(getwd(),fname))
|
|
file.rename(from = file.path(getwd(),fname), to = file.path(pairDirK,fname))
|
|
|
|
#3
|
|
x_rem2_gene <- X[X$NumGenes_Avg_X1 >= 2 & X$NumGenes_Avg_X2 >= 2,]
|
|
plotly_path <- paste(getwd(),"/",outputpath,"/","Scatter_lm_GTF_Averages_",Name1,"_vs_",Name2,"_All_byOverlap_above2genes.html",sep="")
|
|
gg <- ggplot(data = x_rem2_gene,aes(x=Z_lm_K_Avg_X1,y=Z_lm_K_Avg_X2,color=Overlap_Avg,Term=Term_Avg,Genes=Genes_Avg_X1,NumGenes=NumGenes_Avg_X1,AllPossibleGenes=AllPossibleGenes_Avg_X1,SD_1=Z_lm_K_SD_X1,SD_2=Z_lm_K_SD_X2)) +
|
|
xlab(paste("GO Term Avg lm Z for ",Name1,sep="")) +
|
|
geom_rect(aes(xmin=-2,xmax=2,ymin=-2,ymax=2),color="grey20",size=0.25,alpha=0.1,inherit.aes = FALSE,fill=NA) + geom_point(shape=3) +
|
|
ylab(paste("GO Term Avg lm Z for ",Name2,sep="")) + ggtitle(paste("Comparing Average GO Term Z lm for ",Name1," vs. ",Name2,sep="")) +
|
|
theme_Publication_legend_right()
|
|
pdf(paste(getwd(),"/",outputpath,"/","Scatter_lm_GTF_Averages_",Name1,"_vs_",Name2,"_All_ByOverlap_above2genes.pdf",sep=""),width = 12,height = 8)
|
|
gg
|
|
dev.off()
|
|
pgg <- ggplotly(gg)
|
|
#pgg
|
|
print("#3")
|
|
fname <- paste("Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_byOverlap_above2genes.html",sep="")
|
|
print(fname)
|
|
htmlwidgets::saveWidget(pgg, file.path(getwd(),fname))
|
|
file.rename(from = file.path(getwd(),fname), to = file.path(pairDirK,fname))
|
|
|
|
#4
|
|
X_overlap_nothresold <- X[!(is.na(X$Overlap_Avg)),]
|
|
gg <- ggplot(data = X_overlap_nothresold,aes(x=Z_lm_K_Avg_X1,y=Z_lm_K_Avg_X2,color=Overlap_Avg,Term=Term_Avg,Genes=Genes_Avg_X1,NumGenes=NumGenes_Avg_X1,AllPossibleGenes=AllPossibleGenes_Avg_X1,SD_1=Z_lm_K_SD_X1,SD_2=Z_lm_K_SD_X2)) +
|
|
xlab(paste("GO Term Avg lm Z for ",Name1,sep="")) +
|
|
geom_rect(aes(xmin=-2,xmax=2,ymin=-2,ymax=2),color="grey20",size=0.25,alpha=0.1,inherit.aes = FALSE,fill=NA) + geom_point(shape=3) +
|
|
ylab(paste("GO Term Avg lm Z for ",Name2,sep="")) + ggtitle(paste("Comparing Average GO Term Z lm for ",Name1," vs. ",Name2,sep = "")) +
|
|
theme_Publication_legend_right()
|
|
pdf(paste(getwd(),"/",outputpath,"/","Scatter_lm_GTF_Averages_",Name1,"_vs_",Name2,"_Above2SD_ByOverlap.pdf",sep=""),width = 12,height = 8)
|
|
gg
|
|
dev.off()
|
|
pgg <- ggplotly(gg)
|
|
#pgg
|
|
print("#4")
|
|
fname <- paste("Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_Above2SD_ByOverlap.html",sep="")
|
|
print(fname)
|
|
htmlwidgets::saveWidget(pgg, file.path(getwd(),fname))
|
|
file.rename(from = file.path(getwd(),fname), to = file.path(pairDirK,fname))
|
|
|
|
#5
|
|
#only output GTA terms where average score is still above 2 after subtracting the SD
|
|
#Z1 will ID aggravators, Z2 alleviators
|
|
Z1 <- X
|
|
Z1$L_Subtract_SD_X1 <- Z1$Z_lm_K_Avg_X1 - Z1$Z_lm_K_SD_X1
|
|
Z1$L_Subtract_SD_X2 <- Z1$Z_lm_K_Avg_X2 - Z1$Z_lm_K_SD_X2
|
|
|
|
Z2 <- X
|
|
Z2$L_Subtract_SD_X1 <- Z1$Z_lm_K_Avg_X1 + Z1$Z_lm_K_SD_X1
|
|
Z2$L_Subtract_SD_X2 <- Z1$Z_lm_K_Avg_X2 + Z1$Z_lm_K_SD_X2
|
|
|
|
|
|
X1_Specific_Aggravators2 <- Z1[which(Z1$L_Subtract_SD_X1 >= 2 & Z1$L_Subtract_SD_X2 < 2),]
|
|
X1_Specific_Alleviators2 <- Z2[which(Z2$L_Subtract_SD_X1 <= -2 & Z2$L_Subtract_SD_X2 > -2),]
|
|
|
|
X2_Specific_Aggravators2 <- Z1[which(Z1$L_Subtract_SD_X2 >= 2 & Z1$L_Subtract_SD_X1 < 2),]
|
|
X2_Specific_Alleviators2 <- Z2[which(Z2$L_Subtract_SD_X2 <= -2 & Z2$L_Subtract_SD_X1 > -2),]
|
|
|
|
Overlap_Aggravators2 <- Z1[which(Z1$L_Subtract_SD_X1 >= 2 & Z1$L_Subtract_SD_X2 >= 2),]
|
|
Overlap_Alleviators2 <- Z2[which(Z2$L_Subtract_SD_X2 <= -2 & Z2$L_Subtract_SD_X1 <= -2),]
|
|
|
|
X2_Specific_Aggravators2_X1_Alleviatiors2 <- Z1[which(Z1$L_Subtract_SD_X2 >= 2 & Z2$L_Subtract_SD_X1 <= -2),]
|
|
X2_Specific_Alleviators2_X1_Aggravators2 <- Z2[which(Z2$L_Subtract_SD_X2 <= -2 & Z1$L_Subtract_SD_X1 >= 2),]
|
|
|
|
X$Overlap <- NA
|
|
|
|
try(X[X$Term_Avg %in% X1_Specific_Aggravators2$Term_Avg,]$Overlap <- paste(Name1,"Specific_Deletion_Suppressors",sep="_"))
|
|
try(X[X$Term_Avg %in% X1_Specific_Alleviators2$Term_Avg,]$Overlap <- paste(Name1,"Specific_Deletion_Enhancers",sep="_"))
|
|
try(X[X$Term_Avg %in% X2_Specific_Aggravators2$Term_Avg,]$Overlap <- paste(Name2,"Specific_Deletion_Suppressors",sep="_"))
|
|
try(X[X$Term_Avg %in% X2_Specific_Alleviators2$Term_Avg,]$Overlap <- paste(Name2,"Specific_Deletion_Enhancers",sep="_"))
|
|
try(X[X$Term_Avg %in% Overlap_Aggravators2$Term_Avg,]$Overlap <- "Overlapping_Deletion_Suppressors")
|
|
try(X[X$Term_Avg %in% Overlap_Alleviators2$Term_Avg,]$Overlap <- "Overlapping_Deletion_Enhancers")
|
|
try(X[X$Term_Avg %in% X2_Specific_Aggravators2_X1_Alleviatiors2$Term_Avg,]$Overlap <- paste(Name2,"Deletion_Suppressors",Name1,"Deletion_Enhancers",sep="_"))
|
|
try(X[X$Term_Avg %in% X2_Specific_Alleviators2_X1_Aggravators2$Term_Avg,]$Overlap <- paste(Name2,"Deletion_Enhancers",Name1,"Deletion_Suppressors",sep="_"))
|
|
|
|
X_abovethreshold <- X[!(is.na(X$Overlap)),]
|
|
gg <- ggplot(data = X_abovethreshold,aes(x=Z_lm_K_Avg_X1,y=Z_lm_K_Avg_X2,color=Overlap,Term=Term_Avg,Genes=Genes_Avg_X1,NumGenes=NumGenes_Avg_X1,AllPossibleGenes=AllPossibleGenes_Avg_X1,SD_1=Z_lm_K_SD_X1,SD_2=Z_lm_K_SD_X2)) +
|
|
xlab(paste("GO Term Avg lm Z for ",Name1,sep="")) +
|
|
geom_rect(aes(xmin=-2,xmax=2,ymin=-2,ymax=2),color="grey20",size=0.25,alpha=0.1,inherit.aes = FALSE,fill=NA) + geom_point(shape=3) +
|
|
ylab(paste("GO Term Avg lm Z for ",Name2,sep="")) + ggtitle(paste("Comparing Average GO Term Z lm for ",Name1," vs. ",Name2,sep="")) +
|
|
theme_Publication_legend_right()
|
|
pdf(paste(getwd(),"/",outputpath,"/","Scatter_lm_GTF_Averages_",Name1,"_vs_",Name2,"_All_ByOverlap_AboveThreshold.pdf",sep=""),width = 12,height = 8)
|
|
gg
|
|
dev.off()
|
|
pgg <- ggplotly(gg)
|
|
#pgg
|
|
print("#5")
|
|
fname <- paste("Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_ByOverlap_AboveThreshold.html",sep="")
|
|
print(fname)
|
|
htmlwidgets::saveWidget(pgg, file.path(getwd(),fname))
|
|
file.rename(from = file.path(getwd(),fname), to = file.path(pairDirK,fname))
|
|
|
|
#6
|
|
gg <- ggplot(data = X_abovethreshold,aes(x=Z_lm_K_Avg_X1,y=Z_lm_K_Avg_X2,color=Overlap,Term=Term_Avg,Genes=Genes_Avg_X1,NumGenes=NumGenes_Avg_X1,AllPossibleGenes=AllPossibleGenes_Avg_X1,SD_1=Z_lm_K_SD_X1,SD_2=Z_lm_K_SD_X2)) +
|
|
xlab(paste("GO Term Avg lm Z for ",Name1,sep="")) +
|
|
geom_text(aes(label=Term_Avg),nudge_y = 0.25,size=2) +
|
|
geom_rect(aes(xmin=-2,xmax=2,ymin=-2,ymax=2),color="grey20",size=0.25,alpha=0.1,inherit.aes = FALSE,fill=NA) + geom_point(shape=3,size=3) +
|
|
ylab(paste("GO Term Avg lm Z for ",Name2,sep="")) + ggtitle(paste("Comparing Average GO Term Z lm for ",Name1," vs. ",Name2,sep="")) +
|
|
theme_Publication_legend_right()
|
|
pdf(paste(getwd(),"/",outputpath,"/","Scatter_lm_GTF_Averages_",Name1,"_vs_",Name2,"_All_ByOverlap_AboveThreshold_names.pdf",sep=""),width = 20,height = 20)
|
|
gg
|
|
dev.off()
|
|
pgg <- ggplotly(gg)
|
|
#pgg
|
|
print("#6")
|
|
fname <- paste("Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_ByOverlap_AboveThreshold_names.html",sep="")
|
|
print(fname)
|
|
htmlwidgets::saveWidget(pgg, file.path(getwd(),fname))
|
|
file.rename(from = file.path(getwd(),fname), to = file.path(pairDirK,fname))
|
|
|
|
#7
|
|
X_abovethreshold$X1_Rank <- NA
|
|
X_abovethreshold$X1_Rank <- rank(-X_abovethreshold$Z_lm_K_Avg_X1,ties.method = "random")
|
|
X_abovethreshold$X2_Rank <- NA
|
|
X_abovethreshold$X2_Rank <- rank(-X_abovethreshold$Z_lm_K_Avg_X2,ties.method = "random")
|
|
|
|
|
|
gg <- ggplot(data = X_abovethreshold,aes(x=Z_lm_K_Avg_X1,y=Z_lm_K_Avg_X2,color=Overlap,Term=Term_Avg,Genes=Genes_Avg_X1,NumGenes=NumGenes_Avg_X1,AllPossibleGenes=AllPossibleGenes_Avg_X1,SD_1=Z_lm_K_SD_X1,SD_2=Z_lm_K_SD_X2)) +
|
|
xlab(paste("GO Term Avg lm Z for ",Name1,sep="")) +
|
|
geom_text(aes(label=X1_Rank),nudge_y = 0.25,size=4) +
|
|
geom_rect(aes(xmin=-2,xmax=2,ymin=-2,ymax=2),color="grey20",size=0.25,alpha=0.1,inherit.aes = FALSE,fill=NA) + geom_point(shape=3,size=3) +
|
|
ylab(paste("GO Term Avg lm Z for ",Name2,sep="")) + ggtitle(paste("Comparing Average GO Term Z lm for ",Name1," vs. ",Name2,sep="")) +
|
|
theme_Publication_legend_right()
|
|
pdf(paste(getwd(),"/",outputpath,"/","Scatter_lm_GTF_Averages_",Name1,"_vs_",Name2,"_All_ByOverlap_AboveThreshold_numberedX1.pdf",sep=""),width = 15,height = 15)
|
|
gg
|
|
dev.off()
|
|
pgg <- ggplotly(gg)
|
|
#pgg
|
|
print("#7")
|
|
fname <- paste("Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_ByOverlap_AboveThreshold_numberedX1.html",sep="")
|
|
print(fname)
|
|
htmlwidgets::saveWidget(pgg, file.path(getwd(),fname))
|
|
file.rename(from = file.path(getwd(),fname), to = file.path(pairDirK,fname))
|
|
|
|
#8
|
|
gg <- ggplot(data = X_abovethreshold,aes(x=Z_lm_K_Avg_X1,y=Z_lm_K_Avg_X2,color=Overlap,Term=Term_Avg,Genes=Genes_Avg_X1,NumGenes=NumGenes_Avg_X1,AllPossibleGenes=AllPossibleGenes_Avg_X1,SD_1=Z_lm_K_SD_X1,SD_2=Z_lm_K_SD_X2)) +
|
|
xlab(paste("GO Term Avg lm Z for ",Name1,sep="")) +
|
|
geom_text(aes(label=X2_Rank),nudge_y = 0.25,size=4) +
|
|
geom_rect(aes(xmin=-2,xmax=2,ymin=-2,ymax=2),color="grey20",size=0.25,alpha=0.1,inherit.aes = FALSE,fill=NA) + geom_point(shape=3,size=3) +
|
|
ylab(paste("GO Term Avg lm Z for ",Name2,sep="")) + ggtitle(paste("Comparing Average GO Term Z lm for ",Name1," vs. ",Name2,sep="")) +
|
|
theme_Publication_legend_right()
|
|
pdf(paste(getwd(),"/",outputpath,"/","Scatter_lm_GTF_Averages_",Name1,"_vs_",Name2,"_All_ByOverlap_AboveThreshold_numberedX2.pdf",sep=""),width = 15,height = 15)
|
|
gg
|
|
dev.off()
|
|
pgg <- ggplotly(gg)
|
|
#pgg
|
|
print("#8")
|
|
fname <- paste("Scatter_lm_GTA_Averages_",Name1,"_vs_",Name2,"_All_ByOverlap_AboveThreshold_numberedX2.html",sep="")
|
|
print(fname)
|
|
htmlwidgets::saveWidget(pgg, file.path(getwd(),fname))
|
|
file.rename(from = file.path(getwd(),fname), to = file.path(pairDirK,fname))
|
|
|
|
print("write csv files")
|
|
write.csv(x=X,file = paste(getwd(),"/",outputpath,"/","All_GTF_Avg_Scores_",Name1,"_vs_",Name2,".csv",sep=""),row.names = FALSE)
|
|
write.csv(x=X_abovethreshold,file = paste(getwd(),"/",outputpath,"/","AboveThreshold_GTF_Avg_Scores_",Name1,"_vs_",Name2,".csv",sep=""),row.names = FALSE)
|
|
|
|
#End of GTA Pairwise compare for K values
|