[226bc8]: / R / part02.R

Download this file

131 lines (112 with data), 4.8 kB

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
################################################################################
# Function which:
# - calculates the correlations between drug response profiles
# - produces correlation matrix (either symmetric or 2 combined*)
#
# * for example, for two sets of samples coming from patients with different
# diagnosis
################################################################################
makeCorrHeatmap = function(mt, mt2=NA, colsc, concNo="one",
ctab=BloodCancerMultiOmics2017::conctab,
dtab=BloodCancerMultiOmics2017::drugs) {
# quiets concerns of R CMD check "no visible binding for global variable"
NameX=NULL; NameY=NULL; Measure=NULL; x=NULL; y=NULL; xend=NULL;
yend=NULL
# Function which reorders the clustered matrix based on
# computed singular value decomposition.
callback = function(hc, mat){
sv = svd(t(mat))$v[,1]
dend = reorder(as.dendrogram(hc), wts = sv)
rv <- as.hclust(dend)
rv
}
# cluster the drugs
mt = mt[apply(mt, 1, sd)!=0, ]
mt = cor(t(mt))
hc = hclust(as.dist(1-mt))
hc = callback(hc, mt)
ord = rev(hc$label[hc$order])
# prepare data in LF
mt = meltWholeDF(mt)
mt$NameX = factor(mt$Y, levels=rev(ord))
if(concNo=="one") {
mt$NameY = factor(dtab[mt$X,"name"], levels=dtab[ord,"name"])
} else if(concNo=="all") {
mt$NameY = factor(giveDrugLabel(mt$X, ctab, dtab),
levels=giveDrugLabel(ord, ctab, dtab))
}
# take care of the secoond heatmap if present
if(!is.na(mt2)[1]) {
# obtain the correlation coefficients
mt2 = mt2[apply(mt2, 1, sd)!=0, ]
mt2 = cor(t(mt2))
# order the matrix in the same way as mt is ordered
mt2 = mt2[levels(mt$NameX), levels(mt$NameX)]
# put 0 in lower triangle of matrix, include the diagonal
mt2[upper.tri(mt2, diag=TRUE)] = NA
mt2 = meltWholeDF(mt2)
# remove NA rows
mt2 = mt2[!is.na(mt2$Measure),]
# insert top half of matrix data to plotting data frame
idx = match(paste(mt2$X, mt2$Y, sep="."), paste(mt$X, mt$Y, sep="."))
mt[idx, "Measure"] = mt2$Measure
}
# color the drug names on the axis
pathColor = pathColor
pathColor["Other"] = "black"
if(concNo=="one") {
drcol = unname(pathColor[match(dtab[rev(levels(mt$NameX)),"pathway"],
names(pathColor))])
} else if(concNo=="all") {
drcol = unname(pathColor[match(dtab[rev(substring(levels(mt$NameX), 1, 5)),
"pathway"], names(pathColor))])
}
# construct the main plot
mainPlot = ggplot() +
geom_tile(data=mt, aes(x=NameX, y=NameY, fill=Measure)) +
theme_bw() +
theme(axis.text.y = element_text(size=14, colour=drcol),
panel.background=element_rect(fill=NA),
legend.title=element_text(size=14),
legend.text=element_text(size=14)) +
scale_fill_gradientn(name="Correlation coefficient", limits=c(-1,1),
colours=colsc) +
xlab("") + ylab("") + coord_fixed()
# construct the dendrogram
dhc = as.dendrogram(hc)
ddata = dendro_data(dhc, type = "rectangle")
ddata$segments = ddata$segments - 0.5
nXY = nrow(ddata$labels)
dendPlot = ggplot(segment(ddata)) +
geom_segment(aes(x=x, y=y, xend=xend, yend=yend)) +
theme_dendro() +
scale_x_continuous(limits=c(0, nXY), expand = c(0,0)) +
scale_y_reverse()
# construct annotation for drug labels
annoPlot = ggplot(data.frame(x=1, y=factor(names(pathColor),
levels=names(pathColor))),
aes(x, y, fill=y)) +
geom_tile() +
scale_fill_manual("Targeted pathway", values=pathColor) +
theme(legend.title=element_text(size=14), legend.text=element_text(size=14))
# construct gtable
wdths = c(3, 0.22*nXY, 0.1)
hghts = c(0.1, 0.22*nXY, 1)
gt = gtable(widths=unit(wdths, "in"), heights=unit(hghts, "in"))
# make grobs
heat = ggplotGrob(mainPlot)
denX = ggplotGrob(dendPlot)
anno = ggplotGrob(annoPlot)
# fill in the gtable
gt = gtable_add_grob(gt, gtable_filter(heat, "panel"), 2, 2)
gt = gtable_add_grob(gt, denX$grobs[[whichInGrob(denX, "panel")]], 3, 2)
gt = gtable_add_grob(gt, heat$grobs[[whichInGrob(heat, "axis-l")]], 2, 1)
# gtable with legends
wdthsl = c(3, 3)
hghtsl = c(2.5)
gtl = gtable(widths=unit(wdthsl, "in"), heights=unit(hghtsl, "in"))
gtl = gtable_add_grob(gtl, gtable_filter(heat, "guide-box"), 1, 1)
gtl = gtable_add_grob(gtl, gtable_filter(anno, "guide-box"), 1, 2)
return(list("figure"=list(width=sum(wdths), height=sum(hghts), plot=gt),
"legend"=list(width=sum(wdthsl), height=sum(hghtsl), plot=gtl)))
}