# hetset visualization tools ==================================================
# plot heterogeneous networks -------------------------------------------------
plot_hetset <- function(H,plot_in = "no_file",file_name = "hetset_plot_temp"){
if(is.na(H@metadata$slf[1])){
H@metadata$slf <- H@NAMES
H@metadata$color_type <- "raw"
}
plot_type <- as.character(min(3,sum(!is.na(H@metadata$slf))))
if (sum(!is.na(H@metadata$slf))>10){
plot_type <- "make_table"
}
off_flag <- F
# plot on screen or file (select type)
switch(plot_in,
"png" = {
off_flag <- T
png(filename = paste0(file_name,".png"),width = 13,height = 13,
units = "cm",res = 90)
},
"svg" = {
off_flag <- T
svg(filename = paste0(file_name,".svg"),width = 8,height = 8)
},
"pdf" = {
off_flag <- T
pdf(file = paste0(file_name,".pdf"),width = 8,height = 8)
})
# set colors for plotting
color_sel <- colors()[c(130,35,17,50,500,578,610,506,640,636,654,96,24)]
# format labelling and coloring information
if (is.null(H@metadata$color_type)){ # case: "no type specified"
if (!is.null(H$sample_color)){ # subcase: "labels assigned manually"
H@metadata$color_type <- "indiv"
H$sample_color <- as.factor(H$sample_color)
} else {
if (is.null(H$prt)){ # subcase: "raw data"
H@metadata$color_type <- "raw"
H$sample_color <- NA
} else { # subcase: "labels to subpopulations"
H@metadata$color_type <- "subpop"
H$sample_color <- H$prt
}
}
} else { # case: "type specified"
if (H@metadata$color_type %in% c("subpop","pop","prt","raw")){
if (H@metadata$color_type == "raw"){ # subcase: "raw"
H$sample_color <- rep(NA,ncol(H))
} else { # subcase: "subpop"
H@metadata$color_type <- "subpop"
H$sample_color <- H$prt
}
} else {
if (H@metadata$color_type != "indiv"){
# subcase: name of labeling factor in color_type
label <- H@metadata$color_type
eval(parse(text = c("H$sample_color <- H$",label)))
H$sample_color <- as.factor(H$sample_color)
rm(label)
H@metadata$color_type <- "indiv"
} else {
H$sample_color <- as.factor(H$sample_color)
}
}
}
plot.info <- .prep_legend(H = H,color = color_sel)
# type of plot
switch(plot_type,
"1" = {
.plot_hetset_histogram(plot.info)
},
"2" = {
.plot_hetset_scatterplot(plot.info)
},
"3" = {
.plot_hetset_scattermatrix(plot.info)
},
"make_table" = {
cat(" report for more than 10 selected features \n")
cat(" ... not yet implemented ... \n")
})
if (off_flag){
dev.off()
} else {
par(mfrow = c(1,1),mar = c(5,4,4,2)+0.1,mgp = c(3,1,0))
}
}
# prepare legend, colors and plotting order -----------------------------------
.prep_legend <- function(H,color){
na.group.term <- "unlabelled"
# legend for 1 selected feature
if (sum(!is.na(H@metadata$slf)) == 1){
switch(H@metadata$color_type,
"raw" = {
legend_info <- data.frame("typ" = "H.A",
"ID" = NA,
"tex" = "all samples",
"pch" = NA,
"lwd" = NA,
"col" = NA,
"fil" = color[3],
"lty" = NA,
"bor" = "black")
H$sample_color <- 1
},
"subpop" = {
legend_text <- c("all samples",
"subpopulation A",
"subpopulation B",
na.group.term)
legend_info <- data.frame("typ" = c("H.A","H.S","H.S","H.S"),
"ID" = c(NA,"A","B",na.group.term),
"tex" = legend_text,
"pch" = NA,
"lwd" = NA,
"col" = c("white",color[1:3]),
"fil" = c("white",color[1:3]),
"lty" = NA,
"bor" = "black")
H$sample_color <- factor(x = H$sample_color,
levels = c(levels(x = H$sample_color),
na.group.term))
if (sum(is.na(H$sample_color))==0){
legend_info <- legend_info[-4,]
} else {
H$sample_color[is.na(H$sample_color)] <- na.group.term
}
legend_info$ID <- droplevels(legend_info$ID)
H$sample_color <- droplevels(H$sample_color)
},
"indiv" = {
color <- color[-(1:2)]
gr.distr <- table(H$sample_color)
Gi <- names(-sort(-gr.distr))
H$sample_color <- factor(x = H$sample_color,
levels = c(levels(x = H$sample_color),
na.group.term))
H$sample_color[is.na(H$sample_color)] <- na.group.term
Gi <- c(Gi,na.group.term)
legend_text <- c("all samples",Gi)
legend_info <- data.frame("typ" = c("H.A",rep("H.S",length(Gi))),
"ID" = c(NA,Gi),
"tex" = legend_text,
"pch" = NA,
"lwd" = NA,
"col" = NA,
"fil" = c("white",color[2:length(Gi)],
color[1]),
"lty" = NA,
"bor" = "black")
rm(gr.distr,color,Gi,legend_text)
if (sum(H$sample_color==na.group.term)==0){
legend_info <- legend_info[-nrow(legend_info),]
}
H$sample_color <- factor(H$sample_color,
levels = legend_info$"ID"[2:nrow(legend_info)])
legend_info$ID <- droplevels(legend_info$ID)
})
} else { # legend for scatterplots
no.dens <- is.na(H@metadata$prm.A$mean[1])
switch(H@metadata$color_type,
"raw" = {
legend_text <- c(" ","location subpop. A",
"location subpop. B","location all samples")
legend_info <- data.frame("typ" = c("H.A","d","d","d"),
"ID" = c(na.group.term,NA,NA,NA),
"tex" = legend_text,
"pch" = c(16,NA,NA,NA),
"lwd" = c(NA,2,2,2),
"col" = color[c(3,1,2,3)],
"fil" = NA, #color[c(3,1,2,3)],
"lty" = c(0,1,1,1),
"bor" = NA)
H$sample_color <- as.factor(na.group.term)
if (no.dens){
legend_info <- legend_info[1,]
}
legend_info$ID <- droplevels(legend_info$ID)
},
"subpop" = {
legend_text <- c("subpopulation A",
"subpopulation B",
na.group.term,
"location subpop. A",
"location subpop. B",
"location all samples")
legend_info <- data.frame("typ" = c("H.S","H.S","H.S","d","d","d"),
"ID" = c("A","B",na.group.term,NA,NA,NA),
"tex" = legend_text,
"pch" = c(rep(16,3),NA,NA,NA),
"lwd" = c(rep(NA,3),2,2,2),
"col" = color[c(1,2,3,1,2,3)],
"fil" = NA,
"lty" = c(rep(NA,3),1,1,1),
"bor" = NA)
if (no.dens){
legend_info <- legend_info[1:3,]
}
H$sample_color <- factor(x = H$sample_color,
levels = c(levels(x = H$sample_color),
na.group.term))
H$sample_color[is.na(H$sample_color)] <- na.group.term
if (sum(H$sample_color==na.group.term)==0){
legend_info <- legend_info[-3,]
}
H$sample_color <- factor(H$sample_color,
levels = legend_info$"ID"[1:sum(legend_info$typ != "d")])
legend_info$ID <- droplevels(legend_info$ID)
},
"indiv" = {
color.loc <- color[1:3]
color <- color[-(1:2)]
gr.distr <- table(H$sample_color)
Gi <- names(-sort(-gr.distr))
H$sample_color <- factor(x = H$sample_color,
levels = c(levels(x = H$sample_color),
na.group.term))
H$sample_color[is.na(H$sample_color)] <- na.group.term
Gi <- c(Gi,na.group.term)
legend_text <- c(Gi,"location subpop. A","location subpop. B",
"location all samples")
lGi <- length(Gi)
legend_info <- data.frame("typ" = c(rep("H.S",lGi),"d","d","d"),
"ID" = c(Gi,NA,NA,NA),
"tex" = legend_text,
"pch" = c(rep(16,lGi),NA,NA,NA),
"lwd" = c(rep(NA,lGi),2,2,2),
"col" = c(color[2:lGi],1,color.loc),
"fil" = NA,
"lty" = c(rep(0,lGi),1,1,1),
"bor" = NA)
rm(gr.distr,color,color.loc,Gi,legend_text)
if (no.dens){
legend_info <- legend_info[1:lGi,]
}
if (sum(H$sample_color==na.group.term)==0){
legend_info <- legend_info[-lGi,]
}
H$sample_color <- factor(H$sample_color,
levels = legend_info$"ID"[1:nrow(legend_info)])
legend_info$ID <- droplevels(legend_info$ID)
})
}
return(list("H" = H,"legend_info" = legend_info))
}
# plot histogram --------------------------------------------------------------
.plot_hetset_histogram <- function(plot.info){
H <- plot.info$H
info <- plot.info$legend_info
rm(plot.info)
x <- t(assays(H[H@metadata$slf,])[[1]])
x_min <- min(x,na.rm = T)
x_max <- max(x,na.rm = T)
x_breaks <- seq(x_min,x_max,(x_max-x_min)/20)
hist_out <- hist(x,breaks = x_breaks,
main = paste0("Histogram of ",H@metadata$slf),
xlab = paste0("log(",H@metadata$slf,")"))
legend_fill <- rgb(t(col2rgb("white")),maxColorValue = 255)
if (nrow(info) > 1){
for (i in 2:nrow(info)){
Gh <- H$sample_color==info$ID[i]
hist(x[Gh],breaks = x_breaks,add = T,
col = rgb(t(col2rgb(info$col[i])),alpha = 125,maxColorValue = 255),
ann=F,xaxt='n',yaxt='n')
legend_fill <- c(legend_fill,rgb(t(col2rgb(info$col[i])),alpha = 125,
maxColorValue = 255))
}
}
legend("topright",legend = info$tex,bty = 'n',fill = legend_fill)
}
# make scatterplot ------------------------------------------------------------
.plot_hetset_scatterplot <- function(plot.info){
H <- plot.info$H
info <- plot.info$legend_info
rm(plot.info)
# set layout
layout(matrix(c(2,0,1,3),2,2,byrow = T),widths = c(3,1),heights = c(1,3))
par(mar=c(2,2,0,0)+0.1)
# calculate adjustment-points for plot modules
n_bars <- 20
x_barwidth <- diff(range(assays(H[H@metadata$slf[1],])[[1]],na.rm = T))/n_bars
x_lowadj <- min(assays(H[H@metadata$slf[1],])[[1]],na.rm = T)-x_barwidth/2
x_uppadj <- max(assays(H[H@metadata$slf[1],])[[1]],na.rm = T)+x_barwidth/2
y_barwidth <- diff(range(assays(H[H@metadata$slf[2],])[[1]],na.rm = T))/n_bars
y_lowadj <- min(assays(H[H@metadata$slf[2],])[[1]],na.rm = T)-y_barwidth/2
y_uppadj <- max(assays(H[H@metadata$slf[2],])[[1]],na.rm = T)+y_barwidth/2
# scatterplot
N.cols <- sum(info$typ!="d")
plot(x = assays(H[H@metadata$slf[1],])[[1]],
y = assays(H[H@metadata$slf[2],])[[1]],
type = 'o',xlim = c(x_lowadj,x_uppadj),ylim = c(y_lowadj,y_uppadj),
main = NULL,xlab = 'n',ylab = 'n',ann = F,
pch = info$pch[1:N.cols],lty = 0,
col = rgb(t(col2rgb(info$col[H$sample_color])),alpha = 125,
maxColorValue = 255))
# add densities
if (sum(info$typ=="d")>0){
.add_ellipses(He = H,info$col[info$"typ"=="d"])
}
# histograms x- and y-data --------------------------------------------------
.add_histogram(H,var_no = 1,adj = c(x_lowadj,x_uppadj),n_bars = n_bars,
info = info)
.add_histogram(H,var_no = 2,adj = c(y_lowadj,y_uppadj),n_bars = n_bars,
info = info)
}
# add histograms to hetset-scatterplot ----------------------------------------
.add_histogram <- function(H,var_no,adj,n_bars,info){
if(var_no==1){
par(mar=c(0,2,0.42,0)+0.1)
} else {
par(mar=c(2,0,0,0.42)+0.1)
}
m_breaks <- seq(adj[1],adj[2],diff(adj)/n_bars)
hist_out <- hist(assays(H[H@metadata$slf[var_no],])[[1]],
breaks = m_breaks,plot=F)
barplot(height = hist_out$counts,horiz = (var_no==2),space = 0,col="white",
axes = F)
title(main = H@metadata$slf[var_no],line = -1)
for (i in 1:sum(info$typ=="H.S")){
Gh <- H$sample_color==info$ID[i]
barplot(height = hist(assays(H[H@metadata$slf[var_no],Gh])[[1]],
breaks = m_breaks,plot = F)$counts,
horiz = (var_no==2),space = 0,axes = F,
col=rgb(t(col2rgb(info$col[i])),alpha = 125,maxColorValue = 255),
add=T,ann=F,xaxt='n',yaxt='n')
}
if(var_no==1){
legend("topright",legend = info$tex,box.lwd = 0,
bg = rgb(t(col2rgb("white")),alpha = 128,maxColorValue = 255),
col = rgb(t(col2rgb(info$col)),alpha = 175,maxColorValue = 255),
border = info$bor,
pch = info$pch,lty = info$lty)
}
}
# make matrix of scatterplots -------------------------------------------------
.plot_hetset_scattermatrix <- function(plot.info){
H <- plot.info$H
info <- plot.info$legend_info
rm(plot.info)
N.cols <- sum(info$typ!="d")
# organize plot
paare <- t(combn(x = H@metadata$slf,m = 2))
n_paare <- dim(paare)[1]
par(mfrow = c(rep(ceiling(sqrt(n_paare+1)),2)))
par(mgp = c(-1.2,0,0))
par(mar = rep(0.2,4))
for(i in 1:n_paare){
namen <- paare[i,]
x_lowadj <- min(assays(H[namen[1],])[[1]],na.rm = T)
x_uppadj <- max(assays(H[namen[1],])[[1]],na.rm = T)
y_lowadj <- min(assays(H[namen[2],])[[1]],na.rm = T)
y_uppadj <- max(assays(H[namen[2],])[[1]],na.rm = T)
plot(x = assays(H[namen[1],])[[1]],
y = assays(H[namen[2],])[[1]],
type = 'o',xlim = c(x_lowadj,x_uppadj),ylim = c(y_lowadj,y_uppadj),
main = NULL,xlab = namen[1],ylab = namen[2],ann = T,cex.lab = 1.2,
xaxt = 'n',yaxt = 'n',pch = info$pch[1:N.cols],lty = 0,
col = rgb(t(col2rgb(info$col[H$sample_color])),alpha = 125,
maxColorValue = 255))
par(new=F)
# add densities
if (sum(info$typ=="d")>0){
.add_ellipses(He = subset_hetset(H = H,keep.features = namen),
color = info$col[info$"typ"=="d"])
}
}
plot(x = c(0,1),y = c(0,1),xaxt='n',yaxt='n',type = 'n',ann=F,axes=F)
legend("topright",legend = info$tex,bty = 'n',
col = rgb(t(col2rgb(info$col)),alpha = 175,maxColorValue = 255),
border = info$bor,
fill = as.character(info$fil),
lty = info$lty,
pch = info$pch)
}
# add ellipses and correlation values to scatterplots -------------------------
.add_ellipses <- function(He,color){
subpops <- c("A","B")
for (i in 1:length(subpops)){
X <- assays(He[He@metadata$slf[1],He$prt==subpops[i]])[[1]]
mX <- mean(X)
Y <- assays(He[He@metadata$slf[2],He$prt==subpops[i]])[[1]]
mY <- mean(Y)
if ((sd(X,na.rm = T) > 0) & (sd(Y,na.rm = T) > 0)){
cor_out <- cor.test(x = X,y = Y)
if (cor_out$p.value<0.05){
text(x = mX,y = mY,
labels = round(cor_out$estimate,digits = 2),col = "black")
}
mixtools::ellipse(mu = c(mX,mY),
sigma = cov(x = cbind(t(X),t(Y)),use = "complete.obs"),
alpha = 0.5,newplot = F,draw = T,
col = as.character(color[i]),lwd=2)
}
}
if ((sd(assays(He[He@metadata$slf[1],])[[1]]) > 0) &
(sd(assays(He[He@metadata$slf[2],])[[1]]) > 0)){
corF_out <- cor.test(x = assays(He[He@metadata$slf[1],])[[1]],
y = assays(He[He@metadata$slf[2],])[[1]])
if (corF_out$p.value<0.05){
text(x = He@metadata$prm.full$mean[1],y = He@metadata$prm.full$mean[2],
labels = round(corF_out$estimate,digits = 2),col = "black")
}
mixtools::ellipse(mu = He@metadata$prm.full$mean,
sigma = He@metadata$prm.full$cov,
alpha = 0.5,newplot = F,draw = T,
col = as.character(color[3]),lty=1,lwd=1)
}
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.