new supp figures added - wgcna comparisons, fnr

leekgroup · Feb 12, 2019 · 6ca7d93 · 6ca7d93
1 parent 63f3dbf
commit 6ca7d93
Show file tree

Hide file tree

Showing 14 changed files with 594 additions and 0 deletions.
diff --git a/publication_figures/src/makefig3old.R b/publication_figures/src/makefig3old.R
@@ -0,0 +1,68 @@
+library(cowplot)
+library(ggplot2)
+theme_set(theme_cowplot(font_size=9))
+lambda <- seq(0.3,1.0, length.out = 50)
+cat.plot <- c("uncorrected", "RIN", "multi-covariate", "PC")
+# select category to plot
+select_category <- function(category_name, pr_table){
+  pr_table <- pr_table[which(pr_table$type %in% category_name),]
+  pr_table$type <- factor(pr_table$type, levels = category_name)
+  pr_table$lambda <- rep(lambda,length(category_name))
+  pr_table
+}
+plot.thyroid <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_glasso_canonical_thyroid.Rds")
+plot.thyroid <- select_category(cat.plot, plot.thyroid)
+plot.thyroid <- ggplot(plot.thyroid, aes(x = lambda, y = 1-precision, colour = type)) + geom_point(size =0.3) +
+  xlab("lambda") + ylab("FDR")+ggtitle("Thyroid") 
+plot.lung <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_glasso_canonical_lung.Rds")
+plot.lung <- select_category(cat.plot, plot.lung)
+plot.lung <- ggplot(plot.lung, aes(x = lambda, y = 1-precision, colour = type)) + geom_point(size = 0.3) +
+  xlab("lambda") + ylab("FDR")+ggtitle("Lung")
+plot.blood <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_glasso_canonical_blood.Rds")
+plot.blood <- select_category(cat.plot, plot.blood)
+plot.blood <- ggplot(plot.blood, aes(x = lambda, y = 1-precision, colour = type)) + geom_point(size = 0.3) +
+  xlab("lambda") + ylab("FDR")+ggtitle("Whole Blood")
+plot.skin <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_glasso_canonical_skin.Rds")
+plot.skin <- select_category(cat.plot, plot.skin)
+plot.skin <- ggplot(plot.skin, aes(x = lambda, y = 1-precision, colour = type)) + geom_point(size = 0.3) +
+  xlab("lambda") + ylab("FDR")+ggtitle("skin")
+plot.muscle <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_glasso_canonical_muscle.Rds")
+plot.muscle <- select_category(cat.plot, plot.muscle)
+plot.muscle <- ggplot(plot.muscle, aes(x = lambda, y = 1-precision, colour = type)) + geom_point(size = 0.3) +
+  xlab("lambda") + ylab("FDR")+ggtitle("muscle")
+plot.subcutaneous <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_glasso_canonical_subcutaneous.Rds")
+plot.subcutaneous <- select_category(cat.plot, plot.subcutaneous)
+plot.subcutaneous <- ggplot(plot.subcutaneous, aes(x = lambda, y = 1-precision, colour = type)) + geom_point(size = 0.3) +
+  xlab("lambda") + ylab("FDR")+ggtitle("subcutaneous")
+plot.artery <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_glasso_canonical_artery.Rds")
+plot.artery <- select_category(cat.plot, plot.artery)
+plot.artery <- ggplot(plot.artery, aes(x = lambda, y = 1-precision, colour = type)) + geom_point(size = 0.3) +
+  xlab("lambda") + ylab("FDR")+ggtitle("artery")
+plot.nerve <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_glasso_canonical_nerve.Rds")
+plot.nerve <- select_category(cat.plot, plot.nerve)
+plot.nerve <- ggplot(plot.nerve, aes(x = lambda, y = 1-precision, colour = type)) + geom_point(size = 0.3) +
+  xlab("lambda") + ylab("FDR")+ggtitle("nerve")
+###############################
+fig3 <- plot_grid(plot.blood + xlim(0.3,1.0) + ylim(0, 1) + theme(legend.position="none"),
+                  plot.lung + xlim(0.3,1.0) + ylim(0, 1) + theme(legend.position="none"),
+                  # plot.muscle + xlim(0.3,1.0) + ylim(0, 1) + theme(legend.position="none"),
+                  # plot.artery + xlim(0.3,1.0) + ylim(0, 1) + theme(legend.position="none"),
+                  # plot.skin + xlim(0.3,1.0) + ylim(0, 1) + theme(legend.position="none"),
+                  # plot.nerve + xlim(0.3,1.0) + ylim(0, 1) + theme(legend.position="none"),
+                  # plot.subcutaneous + xlim(0.3,1.0) + ylim(0, 1) + theme(legend.position="none"),
+                  plot.thyroid + xlim(0.3,1.0) + ylim(0, 1) +  theme(legend.position="none"),
+                  align = 'vh',
+                  labels = c("a", "b", "c"),
+                  hjust = -1,
+                  nrow = 1
+)
+legend <- get_legend(plot.lung +
+                       theme(legend.key = element_rect(color = "black", linetype = "solid", size = 0.5),
+                             legend.key.size = unit(0.3, "cm"), legend.key.height=unit(1.5,"line")) +
+                       guides(colour = guide_legend(override.aes = list(size= 1))))
+fig3a <- plot_grid( fig3, legend, rel_widths = c(3, .4))
+
+pdf("/work-zfs/abattle4/parsana/networks_correction/publication_figures/fig3.pdf", height = 2.5, width = 7.2)
+print(fig3a)
+dev.off()
+
diff --git a/publication_figures/src/supp_fig13.R b/publication_figures/src/supp_fig13.R
@@ -0,0 +1,86 @@
+library(cowplot)
+library(ggplot2)
+library(parallel)
+library(igraph)
+library(reshape2)
+theme_set(theme_cowplot(font_size=9))
+
+create.igraph <- function(network){
+    # network[which(network != 0)] <- 1
+    diag(network) <- 0
+    network[network!=0] <- 1
+    graph_from_adjacency_matrix(network, mode = "undirected")
+}
+
+
+fn <- dir("../../networks/", recursive=T, full.names = T)
+fn <- fn[grep("glasso", fn)]
+fn <- fn[grep(paste("raw", "rin", "mc", "/pc/", sep = "|"), fn)]
+
+igraph_net <- mclapply(fn, function(x){
+        print(x)
+        load(x)
+        net.igraph <- lapply(dat.net[1:10], function(eachNet){
+                eachNet[lower.tri(eachNet)] <- t(eachNet)[lower.tri(eachNet)]
+                create.igraph(eachNet)
+                })
+        net.igraph
+        }, mc.cores = 16)
+
+# igraph_net <- lapply(igraph_net, function(x,y){
+#         y <- y
+#         this_net <- lapply(x, function(thisx, thisy){
+#                 V(thisx)$name <- thisy
+#                 thisx
+#                 }, y)
+#         this_net
+#         }, dat.gene.symbol$gene_symbol)
+
+
+nodes_with_20 <- mclapply(igraph_net, function(x)
+                sapply(x, function(y) sum(degree(y) >=20)), mc.cores = 16)
+
+names(nodes_with_20) <- fn
+
+raw <- nodes_with_20[grep("raw", fn)]
+rin <- nodes_with_20[grep("rin", fn)]
+mc <- nodes_with_20[grep("mc", fn)]
+pc <- nodes_with_20[grep("/pc/", fn)]
+
+hub_net <- melt(data.frame(raw, rin, mc, pc))
+
+hub_net$tissue <- NA
+hub_net$tissue[grep("Artery", hub_net$variable)] <- "Tibial Artery"
+hub_net$tissue[grep("Blood", hub_net$variable)] <- "Whole Blood"
+hub_net$tissue[grep("Lung", hub_net$variable)] <- "Lung"
+hub_net$tissue[grep("Muscle", hub_net$variable)] <- "Muscle"
+hub_net$tissue[grep("Subcutaneous", hub_net$variable)] <- "Adipose Subcutaneous"
+hub_net$tissue[grep("Nerve", hub_net$variable)] <- "Nerve Tibial"
+hub_net$tissue[grep("Thyroid", hub_net$variable)] <- "Thyroid"
+hub_net$tissue[grep("Skin", hub_net$variable)] <- "Sun exposed Skin"
+
+hub_net$version <- NA
+hub_net$version[grep("raw", hub_net$variable)] <- "Uncorrected"
+hub_net$version[grep("rin", hub_net$variable)] <- "RIN"
+hub_net$version[grep("mc", hub_net$variable)] <- "multi-covariate"
+hub_net$version[grep("..pc.", hub_net$variable)] <- "PC"
+
+
+select_category <- function(category_name, pr_table){
+  pr_table <- pr_table[which(pr_table$version %in% category_name),]
+  pr_table$version <- factor(pr_table$version, levels = category_name)
+  pr_table
+}
+
+which_version <- c("Uncorrected", "RIN", "multi-covariate", "PC")
+
+plot_this <- select_category(which_version, hub_net)
+
+thisplot <- ggplot(aes(y = value, x = tissue, fill = version),
+        data = plot_this)+
+geom_boxplot()+ ylab("# Nodes with >=20 neighbors")+ xlab("Tissue") + theme_classic()+theme(axis.text.x = element_text(angle = 90, hjust = 1))
+
+pdf("../suppfig13.pdf", height = 7, width = 12)
+print(thisplot)
+dev.off()
+
diff --git a/publication_figures/src/supp_fig14.R b/publication_figures/src/supp_fig14.R
@@ -0,0 +1,136 @@
+library(cowplot)
+library(ggplot2)
+theme_set(theme_cowplot(font_size=9))
+cutheights <- seq(0.9,1.0, length.out = 50)
+
+cat.plot <- c("uncorrected", "RIN", "multi-covariate", "PC")
+
+# select category to plot
+select_category <- function(category_name, pr_table){
+  pr_table <- pr_table[which(pr_table$type %in% category_name),]
+  pr_table$type <- factor(pr_table$type, levels = category_name)
+  pr_table$cutheights <- rep(cutheights,length(category_name))
+  pr_table$fnr <- pr_table$false_negatives/ ((4978 * (4978 - 1))/2 - pr_table$density)
+  pr_table
+}
+
+plot.thyroid <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_wgcna-signed_canonical_thyroid.Rds")
+plot.thyroid <- select_category(cat.plot, plot.thyroid)
+plot.thyroid <- ggplot(plot.thyroid, aes(x = cutheights, y  = fnr, colour = type)) + geom_point(size =0.3) +
+  xlab("cut-heights") + ylab("FNR")+ggtitle("Thyroid")
+plot.lung <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_wgcna-signed_canonical_lung.Rds")
+plot.lung <- select_category(cat.plot, plot.lung)
+plot.lung <- ggplot(plot.lung, aes(x = cutheights, y  = fnr, colour = type)) + geom_point(size = 0.3) +
+  xlab("cut-heights") + ylab("FNR")+ggtitle("Lung")
+plot.blood <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_wgcna-signed_canonical_blood.Rds")
+plot.blood <- select_category(cat.plot, plot.blood)
+plot.blood <- ggplot(plot.blood, aes(x = cutheights, y  = fnr, colour = type)) + geom_point(size = 0.3) +
+  xlab("cut-heights") + ylab("FNR")+ggtitle("Whole Blood")
+plot.skin <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_wgcna-signed_canonical_skin.Rds")
+plot.skin <- select_category(cat.plot, plot.skin)
+plot.skin <- ggplot(plot.skin, aes(x = cutheights, y  = fnr, colour = type)) + geom_point(size = 0.3) +
+  xlab("cut-heights") + ylab("FNR")+ggtitle("skin")
+plot.muscle <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_wgcna-signed_canonical_muscle.Rds")
+plot.muscle <- select_category(cat.plot, plot.muscle)
+plot.muscle <- ggplot(plot.muscle, aes(x = cutheights, y  = fnr, colour = type)) + geom_point(size = 0.3) +
+  xlab("cut-heights") + ylab("FNR")+ggtitle("muscle")
+plot.subcutaneous <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_wgcna-signed_canonical_subcutaneous.Rds")
+plot.subcutaneous <- select_category(cat.plot, plot.subcutaneous)
+plot.subcutaneous <- ggplot(plot.subcutaneous, aes(x = cutheights, y  = fnr, colour = type)) + geom_point(size = 0.3) +
+  xlab("cut-heights") + ylab("FNR")+ggtitle("subcutaneous")
+plot.artery <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_wgcna-signed_canonical_artery.Rds")
+plot.artery <- select_category(cat.plot, plot.artery)
+plot.artery <- ggplot(plot.artery, aes(x = cutheights, y  = fnr, colour = type)) + geom_point(size = 0.3) +
+  xlab("cut-heights") + ylab("FNR")+ggtitle("artery")
+plot.nerve <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_wgcna-signed_canonical_nerve.Rds")
+plot.nerve <- select_category(cat.plot, plot.nerve)
+plot.nerve <- ggplot(plot.nerve, aes(x = cutheights, y  = fnr, colour = type)) + geom_point(size = 0.3) +
+  xlab("cut-heights") + ylab("FNR")+ggtitle("nerve")
+###############################
+fig2 <- plot_grid(plot.blood + xlim(0.9,1.0) + ylim(0.045, 0.050) + theme(legend.position="none"),
+                  plot.lung + xlim(0.9,1.0) + ylim(0.045, 0.050) + theme(legend.position="none"),
+                  # plot.muscle + xlim(0.9,1.0) + + theme(legend.position="none"),
+                  # plot.artery + xlim(0.9,1.0) + + theme(legend.position="none"),
+                  # plot.skin + xlim(0.9,1.0) + + theme(legend.position="none"),
+                  # plot.nerve + xlim(0.9,1.0) + + theme(legend.position="none"),
+                  # plot.subcutaneous + xlim(0.9,1.0) + + theme(legend.position="none"),
+                  plot.thyroid + xlim(0.9,1.0) + ylim(0.045, 0.050) + theme(legend.position="none"),
+                  align = 'vh',
+                  labels = c("a", "b", "c"),
+                  hjust = -1,
+                  nrow = 1
+)
+legend <- get_legend(plot.lung +
+                       theme(legend.key = element_rect(color = "black", linetype = "solid", size = 0.5),
+                             legend.key.size = unit(0.3, "cm"), legend.key.height=unit(1.5,"line")) +
+                       guides(colour = guide_legend(override.aes = list(size= 1))))
+fig2a <- plot_grid( fig2, legend, rel_widths = c(3, .4))
+
+
+####################################### Fig 2b glasso ###################################
+lambda <- seq(0.3,1.0, length.out = 50)
+cat.plot <- c("uncorrected", "RIN", "multi-covariate", "PC")
+# select category to plot
+select_category <- function(category_name, pr_table){
+  pr_table <- pr_table[which(pr_table$type %in% category_name),]
+  pr_table$type <- factor(pr_table$type, levels = category_name)
+  pr_table$lambda <- rep(lambda,length(category_name))
+  pr_table$fnr <- pr_table$false_negatives/ ((4978 * (4978 - 1))/2 - pr_table$density)
+  pr_table
+}
+plot.thyroid <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_glasso_canonical_thyroid.Rds")
+plot.thyroid <- select_category(cat.plot, plot.thyroid)
+plot.thyroid <- ggplot(plot.thyroid, aes(x = lambda, y  = fnr, colour = type)) + geom_point(size =0.3) +
+  xlab("lambda") + ylab("FNR")+ggtitle("Thyroid")
+plot.lung <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_glasso_canonical_lung.Rds")
+plot.lung <- select_category(cat.plot, plot.lung)
+plot.lung <- ggplot(plot.lung, aes(x = lambda, y  = fnr, colour = type)) + geom_point(size = 0.3) +
+  xlab("lambda") + ylab("FNR")+ggtitle("Lung")
+plot.blood <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_glasso_canonical_blood.Rds")
+plot.blood <- select_category(cat.plot, plot.blood)
+plot.blood <- ggplot(plot.blood, aes(x = lambda, y  = fnr, colour = type)) + geom_point(size = 0.3) +
+  xlab("lambda") + ylab("FNR")+ggtitle("Whole Blood")
+plot.skin <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_glasso_canonical_skin.Rds")
+plot.skin <- select_category(cat.plot, plot.skin)
+plot.skin <- ggplot(plot.skin, aes(x = lambda, y  = fnr, colour = type)) + geom_point(size = 0.3) +
+  xlab("lambda") + ylab("FNR")+ggtitle("skin")
+plot.muscle <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_glasso_canonical_muscle.Rds")
+plot.muscle <- select_category(cat.plot, plot.muscle)
+plot.muscle <- ggplot(plot.muscle, aes(x = lambda, y  = fnr, colour = type)) + geom_point(size = 0.3) +
+  xlab("lambda") + ylab("FNR")+ggtitle("muscle")
+plot.subcutaneous <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_glasso_canonical_subcutaneous.Rds")
+plot.subcutaneous <- select_category(cat.plot, plot.subcutaneous)
+plot.subcutaneous <- ggplot(plot.subcutaneous, aes(x = lambda, y  = fnr, colour = type)) + geom_point(size = 0.3) +
+  xlab("lambda") + ylab("FNR")+ggtitle("subcutaneous")
+plot.artery <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_glasso_canonical_artery.Rds")
+plot.artery <- select_category(cat.plot, plot.artery)
+plot.artery <- ggplot(plot.artery, aes(x = lambda, y  = fnr, colour = type)) + geom_point(size = 0.3) +
+  xlab("lambda") + ylab("FNR")+ggtitle("artery")
+plot.nerve <- readRDS("/work-zfs/abattle4/parsana/networks_correction/results/PR/pr_density_glasso_canonical_nerve.Rds")
+plot.nerve <- select_category(cat.plot, plot.nerve)
+plot.nerve <- ggplot(plot.nerve, aes(x = lambda, y  = fnr, colour = type)) + geom_point(size = 0.3) +
+  xlab("lambda") + ylab("FNR")+ggtitle("nerve")
+###############################
+fig3 <- plot_grid(plot.blood + xlim(0.3,1.0) + ylim(0.045, 0.050)+ theme(legend.position="none"),
+                  plot.lung + xlim(0.3,1.0) + ylim(0.045, 0.050)+ theme(legend.position="none"),
+                  # plot.muscle + xlim(0.3,1.0) + theme(legend.position="none"),
+                  # plot.artery + xlim(0.3,1.0) + theme(legend.position="none"),
+                  # plot.skin + xlim(0.3,1.0) + theme(legend.position="none"),
+                  # plot.nerve + xlim(0.3,1.0) + theme(legend.position="none"),
+                  # plot.subcutaneous + xlim(0.3,1.0) + theme(legend.position="none"),
+                  plot.thyroid + xlim(0.3,1.0) + ylim(0.045, 0.050) + theme(legend.position="none"),
+                  align = 'vh',
+                  labels = c("a", "b", "c"),
+                  hjust = -1,
+                  nrow = 1
+)
+
+legend <- get_legend(plot.lung +
+                       theme(legend.key = element_rect(color = "black", linetype = "solid", size = 0.5),
+                             legend.key.size = unit(0.3, "cm"), legend.key.height=unit(1.5,"line")) +
+                       guides(colour = guide_legend(override.aes = list(size= 1))))
+fig2b <- plot_grid( fig3, legend, rel_widths = c(3, .4))
+
+pdf("../suppfig14.pdf", height = 5, width = 7.2)
+plot_grid(fig2a, fig2b, nrow = 2)
+dev.off()