# Patrick Cahan (C) 2017
# patrick.cahan@gmail.com
#' @export
sc_violinClass<-function
(sampTab,
classRes,
sid = "cell_name",
dLevel="cluster",
addRand=0,
threshold=0.20,
ncol =1,
sub_cluster = NA){
rownames(sampTab) = sampTab[,sid]
sids <- rownames(sampTab)
colnames(sampTab)[which(colnames(sampTab) == dLevel)] = "cluster"
dLevel = "cluster"
classRes<-classRes[,sids]
stQ2<-cbind(sampTab[sids,], t(classRes[,sids]))
maxX <-apply(classRes, 1, max)
meaVar <-names(which(maxX>threshold))
test <- melt(stQ2, id.vars = c(sid, dLevel), measure.vars = meaVar)
cnames <- colnames(test)
cnames[which(cnames=='value')] <- "classification_score"
cnames[which(cnames=='variable')] <- "cell_type"
colnames(test) <- cnames
xcol = length(unique(sampTab[,dLevel]))
getPalette <- colorRampPalette(brewer.pal(xcol, "Set2"))
if(!is.na(sub_cluster)){
test = test[test$cluster %in% sub_cluster,]
}
ggplot(test, aes(x = cluster, y = classification_score, fill = cluster)) + ylim(0,1) + geom_violin(scale='width', position='dodge', trim=FALSE) +
facet_wrap(~ cell_type, ncol=ncol) +
# scale_y_continuous(
# expand = c(0, 0),
# name = "Class score",
# breaks = c(0, 0.50, 1.0),
# labels = c("0", "0.50", "1.0"),
# limits = c(0,1)
# ) +
coord_cartesian(clip = "off") +
theme_dviz_hgrid() +
theme(
axis.line.x = element_blank(),
axis.ticks.x = element_blank(),
axis.title.y = element_text(size = 8),
axis.title.x = element_text(size = 8),
axis.text.x = element_text(size = 8, angle = 45, hjust=1),
axis.text.y = element_text(size = 8),
legend.title=element_text(size=10),
legend.text=element_text(size=10),
strip.text.x = element_text(size = 8)
) + geom_boxplot(width=0.1, outlier.size=.25) + scale_fill_manual(values = getPalette(xcol))
}
#' @export
prep_umap_class<-function
(classRes,
sampTab,
nrand,
dLevel,
uCates='all',
sid='sample_name',
topPC=10){
stTmp<-addToST(sampTab, nrand=nrand, sid=sid, dLevels=dLevel)
stTmp<-assign_cate(classRes, stTmp)
colnames(stTmp)[2] <- "group"
if(uCates!='all'){
uCates<-unique(stTmp[,"category"])
}
else{
###uCates <- unique(as.vector(sampTab[,dLevel]))
uCates<-rownames(classRes)
}
cat("PCA\n")
pcRes<-prcomp(t(classRes[uCates,]))
if(topPC>length(uCates)){
topPC <- length(uCates)
}
cat("UMAP\n")
uRes<-umap(pcRes$x[,1:topPC], min_dist=.5)
cat("done")
stTmp<-cbind(stTmp, uRes$layout)
colnames(stTmp)[4]<-"umap.1"
colnames(stTmp)[5]<-"umap.2"
stTmp
}
#' @export
plot_umap<-function
(preRes){
ColorRamp <- colorRampPalette(rev(brewer.pal(n = 12,name = "Paired")))(length(unique(preRes$category)))
ggplot(preRes, aes(x=umap.1, y=umap.2, colour=category) ) + geom_point(pch=19, alpha=3/4, size=1) + theme_bw() + scale_colour_manual(values=ColorRamp)
}
#' @export
addToST<-function(sampTab, nrand, sid="sample_name", dLevels=c("description1")){
grpRand<-rep("rand", nrand)
names(grpRand)<-paste("rand_", 1:nrand, sep='')
rTab<-data.frame(sid=names(grpRand))
for(dlev in dLevels){
rTab<-cbind(rTab, grpRand)
}
colnames(rTab)<-c("sid", dLevels)
qTab<-sampTab[,c(sid, dLevels)]
colnames(qTab)[1]<-"sid"
rbind(qTab, rTab)
}
#' @export
assign_cate <- function (classRes, sampTab, cThresh = 0)
{
topCats <- rownames(classRes)[apply(classRes, 2, which.max)]
sampTab <- cbind(sampTab, category = topCats)
sampTab
}
#' @export
get_cate <- function (classRes, sampTab, dLevel, sid, nrand, cThresh=0, keepRand = FALSE)
{
if(is.data.frame(classRes)){
classRes = as.matrix(classRes)
}
if(nrand == 0){
stTmp = sampTab[,c(sid, dLevel)]
}else{
stTmp <- addToST(sampTab, nrand = nrand, sid = sid, dLevels = dLevel)
}
#stTmp <- assign_cate(classRes, stTmp)
colnames(stTmp)[2] <- "group"
topCat_score=c()
topCats = c()
for(i in 1:ncol(classRes)){
tmp = max(classRes[,i])
topCat_score=c(topCat_score, tmp)
if(tmp < cThresh){
tmp2 = "rand"
topCats = c(topCats, tmp2)
}else{
tmp2 = names(classRes[,i][classRes[,i] == tmp])[1]
topCats = c(topCats, tmp2)
}
}
if(keepRand){
sampTab <- cbind(stTmp, category = topCats, scn_score = topCat_score)
return(sampTab)
}
sampTab <- cbind(sampTab, category = topCats[1:nrow(sampTab)], scn_score = topCat_score[1:nrow(sampTab)])
return(sampTab)
}
#' @export
plot_attr <- function (classRes, sampTab, nrand, dLevel, sid = "sample_name", sub_cluster = NA)
{
if (nrand > 0) {
stTmp <- addToST(sampTab, nrand = nrand, sid = sid, dLevels = dLevel)
}
else {
stTmp = sampTab
}
#stTmp <- addToST(sampTab, nrand = nrand, sid = sid, dLevels = dLevel)
stTmp <- assign_cate(classRes, stTmp)
colnames(stTmp)[which(colnames(stTmp) == dLevel)] <- "group"
getPalette = colorRampPalette(brewer.pal(12, "Paired"))
myPal = getPalette(length(unique(stTmp$category)))
p1 = ggplot(stTmp, aes(x = group, fill = category)) + geom_bar(position = "fill",
width = 0.6) + scale_y_continuous(labels = scales::percent) +
scale_fill_manual(values = myPal) + theme_bw() + coord_flip()
if (all(is.na(sub_cluster))) {
p1
}
else {
ggplot(stTmp[which(stTmp$group %in% sub_cluster), ],
aes(x = group, fill = category)) + geom_bar(position = "fill",
width = 0.6) + scale_y_continuous(labels = scales::percent) +
scale_fill_manual(values = myPal) + theme_bw() +
theme(axis.text.x = element_text(angle = 45, vjust = 1,
hjust = 1))
}
}
#' @export
makeColList<-function(sampTab,dLevel,cellName="cell_name",palName="Blues"){
grps<-as.vector(sampTab[,dLevel])
names(grps)<-as.vector(sampTab[,cellName])
cells<-names(grps)
groupNames<-unique(grps)
xcol <- colorRampPalette(rev(brewer.pal(n = 9,name = palName)))(length(groupNames)+1)[1:length(groupNames)]
names(xcol) <- groupNames
anno_colors <- list(group = xcol)
xx<-data.frame(group=as.factor(grps))
rownames(xx)<-cells
list(ann_col = xx, ann_colors = anno_colors)
}
#' @export
hm_sel_col<-function(
expDat,
genes,
clusCols,# from makeColLis
maxPerGrp=100,
cRow=FALSE,
cCol=FALSE,
limits=c(0,10),
toScale=FALSE,
fontsize_row=4){
allgenes<-rownames(expDat)
missingGenes<-setdiff(genes, allgenes)
if(length(missingGenes)>0){
cat("Missing genes: ", paste0(missingGenes, collapse=","), "\n")
genes<-intersect(genes, allgenes)
}
value<-expDat[genes,]
if(toScale){
value <- t(scale(t(value)))
}
value[value < limits[1]] <- limits[1]
value[value > limits[2]] <- limits[2]
anno_colors <- clusCols[["ann_colors"]]
xx <- clusCols[["ann_col"]]
groupNames<-names(anno_colors$group)
cells<-colnames(expDat)
cells2<-vector()
for(groupName in groupNames){
xi<-which(grps==groupName)
if(length(xi)>maxPerGrp){
tmpCells<-sample(cells[xi], maxPerGrp)
}
else{
tmpCells<-cells[xi]
}
cells2<-append(cells2, tmpCells)
}
value<-value[,cells2]
val_col <- colorRampPalette(rev(brewer.pal(n = 12,name = "Spectral")))(25)
pheatmap(value, cluster_rows = cRow, cluster_cols = cCol, color=val_col,
show_colnames = FALSE, annotation_names_row = FALSE,
annotation_col = xx,
annotation_names_col = FALSE, annotation_colors = anno_colors, fontsize_row=fontsize_row)
}
#' Skyline waterfall
#'
#' Skyline waterfall of classification result
#' @param classRes
#' @param cellType
#' @param sampTab
#' @param dLevel
#' @param yval
#'
#' @return nothing
#'
#' @examples
#' skylineClass(classRes, "blood",sampTab, "timepoint", 0.25)
#'
#' @export
skylineClass<-function(
classRes,
cellType,
sampTab,
dLevel,
yval,
cellIdLab="cell_name",
reorder=TRUE){
cellVals<-as.vector(sampTab[,dLevel])
cellAnns<-unique(cellVals)
# assign cols somehow
xdf<-cbind(sampTab, t(classRes))
newXdf<-data.frame();
# re-order by cellAnn then by classification of cellType
if(reorder){
for(group in cellAnns){
xDat<-xdf[which(xdf[,dLevel]==group),]
xDat<-xDat[ order(xDat[,cellType]),]
newXdf<-rbind(newXdf, xDat)
}
}
else{
newXdf<-xdf
}
xi<-which(colnames(newXdf) == cellType)
colnames(newXdf)[xi] <- "vals"
xi<-which(colnames(newXdf) == dLevel)
colnames(newXdf)[xi] <- "group"
newXdf[,cellIdLab]<-factor(newXdf[,cellIdLab], as.vector(newXdf[,cellIdLab]))
xxi<-which(colnames(newXdf)==cellIdLab)
colnames(newXdf)[xxi]<-"cell_name"
#ggplot(data=newXdf, aes(x=cell_name, y=vals, fill=cellAnns)) + geom_bar(stat="identity") + theme_bw()
ggplot(data=newXdf, aes(x=cell_name, y=vals, fill=group)) + geom_bar(stat="identity", width = 1) + theme_bw() + geom_hline(aes(yintercept=yval), colour="#990000", linetype="dashed") + scale_x_discrete(breaks=NULL) + theme(panel.grid.minor=element_blank(),panel.grid.major=element_blank()) + scale_fill_viridis(discrete=TRUE) + ylim(c(0,1))
#newXdf
}
#' heatmap of the enrichment result
#'
#' Heatmap of tthe enrichment result
#' @param esList returned from ks.extract.more
#' @param threshold threshold 0.05 (holm-corrected pval)
#'
#' @return nothing
#'
#' @examples
#' hm_enr(esList)
#'
#' @export
hm_enr<-function
(esList, # has ES and pval
threshold=0.05,
cRows=FALSE,
cCols=FALSE,
fsr=4
){
cools <- colorRampPalette(rev(brewer.pal(n = 11,name = "RdBu")))(100)
# cools<-colorRampPalette(c("blue", "white", "yellow"))( 100 )
pes <- esList[['ES']]
pp <- esList[['pval']]
pes[which(pp>threshold)]<-0
mmin<-min(pes)
mmax<-max(pes)
mmax<-max(abs(mmin), abs(mmax))
bcol<-"white"
pheatmap(pes,
col=cools,
breaks=seq(from= -1 * mmax, to=mmax, length.out=100),
border_color=bcol,
cluster_rows = cRows,
cluster_cols = cCols,
fontsize_row=fsr)
}
#' re-order cells for plotting
#'
#' re-order cells for plotting
#'
#' @param nvect named vector amed vector of cell-> grp
#' @param newOrder vector of grp names in new order
#'
#' @return named vector new order
#'
#' @export
reorderCellsByGrp<-function(
nvect,# named vector of cell-> grp
newOrder # vector of grp names in new order
){
newAns<-vector()
for(i in 1:length(newOrder)){
nname<-newOrder[i]
aib<-nvect[which(nvect==nname)]
cat(nname," ", length(aib),"\n")
newAns<-append(newAns, aib)
}
newAns
}
#This will spit out a matrix of genepairs by (query samples + avg of training data) ordered by importance in the RF model
#' @title
#' Compile genePairs comparison matrix
#' @description
#' compile genePairs comparison matrix for users to elucidate the biological significance behind classification results
#'
#' @param query_exp the expression matrix of query samples
#' @param training_exp the expression matrix of training samples
#' @param training_st the sample table of training data
#' @param classCol the column name of the sample table that indicates the classes
#' @param sampleCol the column name of the sample table that indicates the sample names. NULL if sample names are indcated in rownames of the sample table
#' @param cnProc the cnProc from the training
#' @param numPairs the number of genes to extract for comparison from most important to least important in the classifier
#'
#' @return gene pair comparison matrix
#'
#' @export
compareGenePairs<-function(query_exp, training_exp, training_st, classCol, sampleCol = NULL, RF_classifier, numPairs = 20, trainingOnly = TRUE) {
importantPairs = RF_classifier$importance
importantPairs = sort(importantPairs[, 1], decreasing = TRUE)
importantPairs = importantPairs[grep("_", names(importantPairs))]
if (numPairs > length(importantPairs)) {
userPairs = importantPairs
}else {
userPairs = importantPairs[1:numPairs]
}
query_pairs = query_transform(expDat = query_exp, genePairs = names(userPairs))
training_pairs = query_transform(expDat = training_exp, genePairs = names(userPairs))
avg_training_pairs = avgGeneCat(expDat = training_pairs, sampTab = training_st, dLevel = classCol, sampID = sampleCol)
if(!trainingOnly){
PairCompareMatrix = makeGeneCompareTab(queryExpTab = query_pairs,
avgGeneTab = avg_training_pairs,
geneSamples = names(userPairs))
PairCompareMatrix
} else{
avg_training_pairs
}
}
#this is visualizing the matrix above
#' @title
#' Gene expression plotting
#' @description
#' Plot gene expressions for visualization and comparison
#'
#' @param expDat comparison expression matrix from \code{\link{makeGeneCompareTab}}
#' @param fontsize_row row font size
#' @param cRows TRUE if cluster the rows
#' @param cCols TRUE if cluster the columns
#'
#' @return a heatmap of gene expressions for the comparsion expression matrix
#' @importFrom RColorBrewer brewer.pal
#' @export
plotGeneComparison<-function(expDat, grps, fontsize_row = 6, cRows = FALSE, cCols = FALSE, toScale = FALSE, ...){
#grps = as.vector(colnames(expDat))
#names(grps) = colnames(expDat)
genes = rownames(expDat)
value<-expDat[genes,] #select the matrix with cgenes
if(toScale){
value <- t(scale(t(value))) #scales the matrix
}
limits=c(0,10)
value[value < limits[1]] <- limits[1] # ensures 0 is the smallest
value[value > limits[2]] <- limits[2] # ensures 10 is the highest
groupNames<-unique(grps) # gather the names of the groups
cells<-names(grps)
cells2<-vector()
for(groupName in groupNames){
xi<-which(grps==groupName) #select samples that are in a certain group
tmpCells<-cells[xi] #if not over the maximum number of samples per group, then use all the samples available
cells2<-append(cells2, tmpCells) # create a vector with all the samples selected for plotting
}
value<-value[,cells2] # select the samples that are going to be used for plotting
xcol <- colorRampPalette(rev(brewer.pal(n = 12,name = "Paired")))(length(groupNames))
names(xcol) <- groupNames
anno_colors <- list(group = xcol)
xx<-data.frame(group=as.factor(grps))
rownames(xx)<-cells
# val_col <- colorRampPalette(rev(RColorBrewer::brewer.pal(n = 12,name = "Spectral")))(25)
val_col<-colorRampPalette(rev(RColorBrewer::brewer.pal(n = 12,name = "Spectral")))(25)
pheatmap(value, color=val_col, cluster_row = cRows, cluster_cols = cCols,
show_colnames = FALSE, annotation_names_row = FALSE,
annotation_col = xx,
annotation_names_col = FALSE,
annotation_colors = anno_colors,
fontsize_row=fontsize_row, ...)
}
#' @title
#' Average Category Gene Expressions
#' @description
#' Averaging the gene expression level for each category
#'
#' @param expDat the normalized expression table from \code{\link{trans_prop}}
#' @param sampTab the sample table
#' @param dLevel the name of the column that contains categories
#' @param sampID the name of the column that contains sample ID
#'
#' @return an average matrix of the gene expressions
#' @export
avgGeneCat<-function(expDat, sampTab, dLevel, sampID = NULL){
if (is.null(sampID) == TRUE) {
sampID = "sampID"
sampTab[, sampID] = rownames(sampTab)
}
returnMatrix = matrix(nrow=nrow(expDat), ncol=0)
rownames(returnMatrix) = rownames(expDat)
for (cat in unique(sampTab[, dLevel])) {
tempExpDat = expDat[, as.vector(sampTab[sampTab[, dLevel] == cat, sampID])]
# calculates the mean
tempMatrix = matrix(apply(tempExpDat, 1, mean), ncol = 1, nrow = length(apply(tempExpDat, 1, mean)))
rownames(tempMatrix) = names(apply(tempExpDat, 1, mean))
colnames(tempMatrix) = paste0(cat, "_Avg")
if(all(rownames(returnMatrix) == rownames(tempMatrix))) {
returnMatrix = cbind(returnMatrix, tempMatrix)
}
}
#return
returnMatrix
}
#' @title
#' Make Gene Comparison Table
#' @description
#' To compile an expression table for comparison
#'
#' @param queryExpTab a matrix of normalized expression query data from \code{\link{trans_prop}}
#' @param avgGeneTab a matrix of averaged expression of training data from \code{\link{avgGeneCat}}.
#' @param querySample a vector or string indicating the query samples
#' @param trainingCat a vector or string indicating the categories of the training data
#' @param geneSamples a vector or string indicating the genes of interest
#'
#' @return a matrix that combines query and training data with genes of interest for comparison
#' @export
makeGeneCompareTab<-function(queryExpTab, avgGeneTab, querySample = NULL, trainingCat=NULL, geneSamples) {
if(is.null(querySample) == TRUE) {
filteredQuery = queryExpTab[geneSamples, ]
querySample = colnames(filteredQuery)
}
else if(all(querySample %in% colnames(queryExpTab)) == FALSE) {
filteredQuery = queryExpTab[geneSamples, ]
querySample = colnames(filteredQuery)
print("Please enter sample names that are in the query table")
}
else {
filteredQuery = queryExpTab[geneSamples, querySample]
}
#
if(is.null(trainingCat) == TRUE) {
filteredGeneTab = avgGeneTab[geneSamples, ]
trainingCat = colnames(filteredGeneTab)
}
else if(all(trainingCat %in% colnames(avgGeneTab)) == FALSE) { #maybe modified it later
filteredGeneTab = avgGeneTab[geneSamples, ]
trainingCat = colnames(filteredGeneTab)
print("Please enter proper category names that are in the filtered gene table")
}
else {
filteredGeneTab = avgGeneTab[geneSamples, trainingCat]
}
if(all(rownames(filteredQuery) == rownames(filteredGeneTab)) == TRUE) {
print("All Good")
}
returnLabel = c(querySample, trainingCat)
returnMatrix = cbind(filteredQuery, filteredGeneTab)
colnames(returnMatrix) = returnLabel
#return
returnMatrix
}
#' @title
#' Make training label
#' @description
#' To make name vector for plotting purpose
#'
#' @param gpTab
#' @param stTrain sample table for training data
#' @param dLevel the annotation for cell type label
#'
#' @return a named vector for plotting purposes
#' @export
findAvgLabel <- function(gpTab, stTrain, dLevel){
sla_avg = colnames(gpTab)[(ncol(gpTab) - length(unique(stTrain[,dLevel])) + 1):ncol(gpTab)]
sla = setNames(sla_avg,sla_avg)
}
#' make tsne from pca
#'
#' make tsne from pca
#'
#' @param expDat expDat
#' @param recRes result of running gpaRecurse
#' @param perplexity (30)
#' @param theta (0.30)
#' @param weighted (TRUE) whether to use PCs from deeper lelves, and to weight them
#'
#' @return tsne matrix
#'
#' @export
pca_to_tsne<-function
(expDat,
recRes, # result of running gpaRecurse
perplexity=30,
theta=0.30,
weighted=TRUE){
tmpMat<-pca_project_all(expDat, recRes)
if(weighted){
weights<-rep(1, ncol(tmpMat))
xi<- grep("L1_",colnames(tmpMat))
weights[xi]<-10
xi2<- grep("L2_",colnames(tmpMat))
weights[xi2]<-.5
xi3<- grep("L3_",colnames(tmpMat))
weights[xi3]<-.25
xi4<- grep("L4_",colnames(tmpMat))
weights[xi4]<-.1
datMat <- tmpMat %*% diag(weights)
}
else{
###datMat<-recRes$results[[1]]$gpRes$pcaRes$pcaRes$x
datMat<-tmpMat
}
tres<-Rtsne(datMat, pca=FALSE, perplexity=perplexity, theta=theta, max_iter=2e3)
xres<-tres$Y
colnames(xres)<-c("TSNE.1", "TSNE.2")
rownames(xres)<-rownames(datMat)
xres
}
project_pca<-function
(expDat,
pcRes,
pcs=FALSE){
if(length(pcs)==1){
if(!pcs){
pcs<-1:ncol(pcRes$rotation)
}
}
scale(t(expDat), pcRes$center, pcRes$scale) %*% pcRes$rotation[,pcs]
}
pca_project_gpa<-function
(expDat,
gpaRecRes,
lName){
myResult<-gpaRecRes$results[[lName]]$gpRes$pcaRes
# get the varGenes
varGenes<-myResult$varGenes
# get the PCs
pcs<-1:ncol(myResult$pcaRes$x)
ans<-project_pca(expDat[varGenes,], myResult$pcaRes, pcs=pcs)
ans
}
pca_project_all<-function
(expDat,
gpaRecRes){
resNames<-names(gpaRecRes$results)
tmpAns<-list()
colCount<-0
allCnames<-vector()
for(resName in resNames){
cat(resName,"\n")
tmpAns[[resName]]<-pca_project_gpa(expDat, gpaRecRes, resName)
cnames<-paste(resName, "_",colnames(tmpAns[[resName]]), sep='')
colnames(tmpAns[[resName]])<-cnames
colCount<-colCount + ncol(tmpAns[[resName]])
allCnames<-append(allCnames, cnames)
}
ans<-matrix(0, nrow=ncol(expDat), ncol=colCount)
rownames(ans)<-colnames(expDat)
stp<-0
for(resName in resNames){
str<-stp+1
stp<-str + ncol(tmpAns[[resName]]) - 1
ans[,str:stp]<-tmpAns[[resName]]
}
colnames(ans)<-allCnames
ans
}
myGrpSort<-function(grps){
grpNames<-unique(grps)
### llevels<- as.numeric(unlist(lapply(strsplit(grpNames, "L"), "[[",2))) NEED TO FIX THIS TO repeat over Ls
xix<-as.numeric(unlist(lapply(strsplit(grpNames, "_G"), "[[",2)))
grpNames[order(xix)]
}
#' heatmap genes and groups
#'
#' heatmap genes and groups
#'
#' @param expDat expDat
#' @param genes genes
#' @param grps vector of cellnames -> grp label
#' @param maxPerGrp 100
#' @param cRow =FALSE,
#' @param cCol =FALSE,
#' @param limits =c(0,10),
#' @param toScale =FALSE,
#' @param fontsize_row =4
#'
#' @return pheatmap
#'
#' @export
hm_gpa_sel<-function(
expDat,
genes,
grps, ## vector of cellnames -> grp label
maxPerGrp=100,
cRow=FALSE,
cCol=FALSE,
limits=c(0,10),
toScale=FALSE,
fontsize_row=4,
reOrderCells=FALSE){
allgenes<-rownames(expDat)
missingGenes<-setdiff(genes, allgenes)
if(length(missingGenes)>0){
cat("Missing genes: ", paste0(missingGenes, collapse=","), "\n")
genes<-intersect(genes, allgenes)
}
value<-expDat[genes,]
if(toScale){
if(class(value)[1]!='matrix'){
value <- t(scale(Matrix::t(value)))
}
else{
value <- t(scale(t(value)))
}
}
value[value < limits[1]] <- limits[1]
value[value > limits[2]] <- limits[2]
groupNames<-unique(grps)
if(reOrderCells){
grps<-grps[order(grps)]
groupNames<-sort(unique(grps))
}
cells<-names(grps)
##
## groupNames<-myGrpSort(grps)
##
cells2<-vector()
for(groupName in groupNames){
xi<-which(grps==groupName)
if(length(xi)>maxPerGrp){
tmpCells<-sample(cells[xi], maxPerGrp)
}
else{
tmpCells<-cells[xi]
}
cells2<-append(cells2, tmpCells)
}
value<-value[,cells2]
xcol <- colorRampPalette(rev(brewer.pal(n = 12,name = "Paired")))(length(groupNames))
names(xcol) <- groupNames
anno_colors <- list(group = xcol)
xx<-data.frame(group=as.factor(grps))
rownames(xx)<-cells
val_col <- colorRampPalette(rev(brewer.pal(n = 12,name = "Spectral")))(25)
#val_col <- colorRampPalette(brewer.pal(n = 12,name = "Spectral"))(100)
pheatmap(value, cluster_rows = cRow, cluster_cols = cCol, color=val_col,
show_colnames = FALSE, annotation_names_row = FALSE,
## annotation_col = annTab,
annotation_col = xx,
annotation_names_col = FALSE, annotation_colors = anno_colors, fontsize_row=fontsize_row)
}
reorderCells<-function
(grpList){
curr_grps<-grpList[[1]]
curr_cells<-names(curr_grps)
for(i in 2:length(grpList)){
prior_grps<-curr_grps
prior_cells<-curr_cells
xList<-grpList[[i]]
prior_names<-sort(unique(prior_grps))
for(j in seq(length(prior_names))){
pname<-prior_names[j]
xi<-which(prior_grps==pname)
tmpCells<-prior_cells[xi]
tmpRes<-sort(xList[tmpCells])
curr_cells[xi]<-names(tmpRes)
curr_grps[xi] <-tmpRes
}
}
names(curr_grps)<-curr_cells
curr_grps
}
#' @export
corplot_sub<-function
(gpaRes,
expDat,
prop=0.1,
pSide=FALSE,
minCount=20,
llevel=1){
###orderedCells<-reorderCells(gpaRes$grp_list)
orderedCells<-reorderCells(gpaRes$grp_list)
ssamp<-subsample_min(orderedCells, prop=prop, minCount=minCount)
##ssamp<-subsample_min(gpaRes$groups, prop=prop, minCount=minCount)
### ssamp<-subsample_min(gpaRes$grp_list[[2]], prop=prop, minCount=minCount)
genes<-getVarFromList(gpaRes, llevel=llevel)
###
### xcor<-cor(expDat[genes,names(ssamp)])
###
xcor<- dist(t(expDat[genes,names(ssamp)]))
##xcor<-as.matrix(gpaRes$results[["L1_G1"]]$gpRes$xdist)[names(ssamp), names(ssamp)]
maxX<-max(xcor)
xcor<-as.matrix((maxX - xcor)/maxX)
llevels<-length(gpaRes$grp_list)
xx<-gpaRes$grp_list[[2]][names(ssamp)]
xx<-data.frame(level_1=as.factor(xx))
cnames<-c("level_1")
if(llevels>2){
for(i in 3:llevels){
cnames<-append(cnames, paste0("level_", i-1))
danss<-gpaRes$grp_list[[i]][names(ssamp)]
xx<-cbind(xx, as.factor(danss))
}
colnames(xx)<-cnames
}
if(pSide){
#topGenes<-gpaRes$groupTree$Get("topGenes")
topGenes<-getTopGenesList(gpaRes$results[["L1_G1"]],7)
### xy<-data.frame(levelX=xx[,ncol(xx)], genes=rep("", nrow(xx)))
xy<-data.frame(levelX=xx[,1], genes=rep("", nrow(xx)))
rLabels<-rep("", nrow(xy))
grpNames<-unique(xy$levelX)
for(grpName in grpNames){
cat("***",grpName,"***\n")
xi<-which(xy$levelX==grpName)
coord<-ceiling( (max(xi)-min(xi)) / 2 ) + min(xi)
cat(grpName," xi:",xi[1], "length: ", length(xi), "coord: ", coord,"\n")
rLabels[coord]<-topGenes[grpName]
}
## ColorRamp <- colorRampPalette(brewer.pal(n = 9,name = "Blues"))(100)
rownames(xy)<-rownames(xx)
pheatmap(xcor,
## color=ColorRamp,
cluster_rows = FALSE,
cluster_cols = FALSE,
show_colnames = FALSE,
show_rownames=TRUE,
annotation_names_row = FALSE,
annotation_col = xx,
# annotation_row = xy,
labels_row=rLabels,
fontsize_row=5)
}
else{
pheatmap(xcor,
cluster_rows = FALSE,
cluster_cols = FALSE,
show_colnames = FALSE,
show_rownames=FALSE,
annotation_names_row = FALSE,
annotation_col = xx)
}
}
#' @export
corplot_triangle<-function
(gpaRes,
expDat,
prop=0.1,
pSide=FALSE,
minCount=20,
llevel=1){
###orderedCells<-reorderCells(gpaRes$grp_list)
orderedCells<-reorderCells(gpaRes$grp_list)
ssamp<-subsample_min(orderedCells, prop=prop, minCount=minCount)
##ssamp<-subsample_min(gpaRes$groups, prop=prop, minCount=minCount)
### ssamp<-subsample_min(gpaRes$grp_list[[2]], prop=prop, minCount=minCount)
genes<-getVarFromList(gpaRes, llevel=llevel)
###
### xcor<-cor(expDat[genes,names(ssamp)])
###
### xcor<- dist(t(expDat[genes,names(ssamp)]))
xcor<-as.matrix(gpaRes$results[["L1_G1"]]$gpRes$xdist)[names(ssamp), names(ssamp)]
maxX<-max(xcor)
xcor<-as.matrix((maxX - xcor)/maxX)
xcor[lower.tri(xcor)]<-0
llevels<-length(gpaRes$grp_list)
xx<-gpaRes$grp_list[[2]][names(ssamp)]
xx<-data.frame(level_1=as.factor(xx))
cnames<-c("level_1")
if(llevels>2){
for(i in 3:llevels){
cnames<-append(cnames, paste0("level_", i-1))
danss<-gpaRes$grp_list[[i]][names(ssamp)]
xx<-cbind(xx, as.factor(danss))
}
colnames(xx)<-cnames
}
if(pSide){
#topGenes<-gpaRes$groupTree$Get("topGenes")
topGenes<-getTopGenesList(gpaRes$results[["L1_G1"]],7)
### xy<-data.frame(levelX=xx[,ncol(xx)], genes=rep("", nrow(xx)))
xy<-data.frame(levelX=xx[,1], genes=rep("", nrow(xx)))
rLabels<-rep("", nrow(xy))
grpNames<-unique(xy$levelX)
for(grpName in grpNames){
cat("***",grpName,"***\n")
xi<-which(xy$levelX==grpName)
coord<-ceiling( (max(xi)-min(xi)) / 2 ) + min(xi)
cat(grpName," xi:",xi[1], "length: ", length(xi), "coord: ", coord,"\n")
rLabels[coord]<-topGenes[grpName]
}
ColorRamp <- colorRampPalette(brewer.pal(n = 9,name = "Blues"))(100)
rownames(xy)<-rownames(xx)
pheatmap(xcor,
color=ColorRamp,
cluster_rows = FALSE,
cluster_cols = FALSE,
show_colnames = FALSE,
show_rownames=TRUE,
annotation_names_row = FALSE,
annotation_col = xx,
# annotation_row = xy,
labels_row=rLabels,
fontsize_row=5)
}
else{
pheatmap(xcor,
cluster_rows = FALSE,
cluster_cols = FALSE,
show_colnames = FALSE,
show_rownames=FALSE,
annotation_names_row = FALSE,
annotation_col = xx)
}
}
getVarFromList<-function(
gpaRes,
llevel=1
){
gresNames<-names(gpaRes$results)
charMatch<-paste0("L",llevel,"_")
gresNames<-gresNames[grep(charMatch, gresNames)]
if(FALSE){
if(llevel==1){
tmpAns<-gpaRes$results[[1]]$gpRes$pcaRes$varGenes
}
else{
tmpAns<-vector()
for(i in seq(length(gpaRes$results))){
tmpAns<-append(tmpAns, gpaRes$results[[i]]$gpRes$pcaRes$varGenes)
}
}
}
tmpAns<-vector()
for(gName in gresNames){
tmpAns<-append(tmpAns, gpaRes$results[[gName]]$gpRes$pcaRes$varGenes)
}
sort(unique(tmpAns))
}
subsample_min<-function
(namedVect,
prop=0.10,
minCount=20){
groups<-unique(namedVect)
grpCounts<-table(namedVect)
newVect<-vector()
for(grp in groups){
cat(grp, ": ", ceiling(prop*grpCounts[grp]),"\n")
numToSamp<-max(minCount, ceiling(prop*grpCounts[grp]))
xi<-which(namedVect==grp)
tmpVect<-sample(namedVect[xi], numToSamp)
newVect<-append(newVect, tmpVect)
}
newVect
}
#' heatmap of the classification result
#'
#' Heatmap of the classification result.
#' @param classMat classMat
#' @param isBig is this a big heatmap? TRUE or FALSE
#' @param cluster_cols cluster_cols
#'
#' @return nothing
#'
#' @examples
#' cn_HmClass(cnRes, isBig=TRUE)
#'
#' @export
sc_hmClass<-function(
classMat,
grps, ## vector of cellnames -> grp label
isBig=FALSE,
maxPerGrp=100,
cRow=FALSE,
cCol=FALSE,
fontsize_row=8,
scale=FALSE,
label_reorder = TRUE
){
cools<-colorRampPalette(c("black", "limegreen", "yellow"))( 100 )
bcol<-'white';
if(isBig){
bcol<-NA;
}
if(label_reorder){
grps <- grps[order(grps)]
cells <- names(grps)
groupNames <- sort(unique(grps))
}else{
cells <- names(grps)
groupNames <- unique(grps)
}
cells2<-vector()
for(groupName in groupNames){
xi<-which(grps==groupName)
if(length(xi)>maxPerGrp){
tmpCells<-sample(cells[xi], maxPerGrp)
}
else{
tmpCells<-cells[xi]
}
cells2<-append(cells2, tmpCells)
}
classMat<-classMat[,cells2]
xcol <- colorRampPalette(rev(brewer.pal(n = 12,name = "Paired")))(length(groupNames))
names(xcol) <- groupNames
anno_colors <- list(group = xcol)
xx<-data.frame(group=as.factor(grps))
rownames(xx)<-cells
if(scale){
mymin<-min(classMat)
mymax<-max(classMat)
}
else{
mymin<-0
mymax<-1
}
pheatmap(classMat, col=cools, breaks=seq(from=mymin, to=mymax, length.out=100), cluster_rows = cRow, cluster_cols = cCol,
show_colnames = FALSE, annotation_names_row = FALSE,
## annotation_col = annTab,
annotation_col = xx,
annotation_names_col = FALSE, annotation_colors = anno_colors, fontsize_row=fontsize_row)
}
#' @export
hmcellsgenes<-function
(expDat,
geneList,
sampTab,# must of a colname 'group' corresponding to names of geneList
n_genes=5,
limits=c(-3,3),
fontsize_row=NULL){
dLevel<-"group"
geneList<-geneList[order(as.numeric(names(geneList)))]
groupNames<-names(geneList)
stX<-data.frame()
for(groupName in groupNames){
stX<-rbind(stX, sampTab[sampTab$group==groupName,])
}
stX<-droplevels(stX)
topx<-unlist(lapply(geneList, length))
topx[which(topx>=n_genes)]<-n_genes
# topx<-topx[order(as.numeric(names(topx)))]
geneList<-geneList[ names(topx) ]
genes<-list()
for(groupName in groupNames){
if(topx[groupName]>0){
genes[[groupName]]<-geneList[[groupName]][1:as.vector(topx[groupName])]
}
}
genes<-unlist(genes)
cat("genes: ", length(genes),"\n")
stX<-stX[order(stX[,dLevel]),]
expDat<-expDat[,rownames(stX)]
value <- t(scale(t(expDat[genes,])))
value[value < limits[1]] <- limits[1]
value[value > limits[2]] <- limits[2]
gene_annotation<-data.frame(group=rep(names(topx), topx))
genes2<-make.unique(genes)
rownames(gene_annotation)<-genes2
rownames(value)<-genes2
cell_annotation<-data.frame(group=stX[,dLevel])
rownames(cell_annotation)<-rownames(stX)
xcol <- colorRampPalette(rev(brewer.pal(n = 12,name = "Paired")))(length(groupNames))
names(xcol) <- groupNames
anno_colors <- list(group = xcol)
pheatmap(value, cluster_rows = FALSE, cluster_cols = FALSE,
show_colnames = F, annotation_names_row = FALSE,
annotation_col = cell_annotation, annotation_row = gene_annotation,
annotation_names_col = FALSE, annotation_colors = anno_colors, fontsize_row=fontsize_row);
list(c1=cell_annotation, g1=gene_annotation)
}
#' @export
hmgenesSimple<-function
(expDat,
cRow=FALSE,
cCol=FALSE,
limits=c(-3,3),
toScale=FALSE,
fontsize_row=NULL,
anc=FALSE){
value<-expDat
if(toScale){
value <- t(scale(t(expDat)))
}
value[value < limits[1]] <- limits[1]
value[value > limits[2]] <- limits[2]
pheatmap(value,
cluster_rows = cRow,
cluster_cols = cCol,
show_colnames = anc, annotation_names_row = FALSE, scale='none',
# clustering_distance_rows='euclidean',
clustering_distance_rows='correlation',
clustering_distance_cols="correlation",
#clustering_distance_cols="euclidean",
clustering_method="average",
annotation_names_col = FALSE, fontsize_row=fontsize_row)
# value;
}
#' @export
hmgenesSimple2<-function
(expDat,
annTab,
cName="group",
cRow=FALSE,
cCol=FALSE,
limits=c(-3,3),
toScale=FALSE,
fontsize_row=NULL){
value<-expDat
if(toScale){
value <- t(scale(t(expDat)))
}
value[value < limits[1]] <- limits[1]
value[value > limits[2]] <- limits[2]
groupNames<-unique(annTab[,cName])
xcol <- colorRampPalette(rev(brewer.pal(n = 12,name = "Paired")))(length(groupNames))
names(xcol) <- groupNames
anno_colors <- list(group = xcol)
if(FALSE){
pheatmap(value,
cluster_rows = cRow,
cluster_cols = cCol,
show_colnames = F, annotation_names_row = FALSE, annotation_names_col = TRUE, scale='none',
# clustering_distance_rows='euclidean',
clustering_distance_rows='correlation',
clustering_distance_cols="correlation",
#clustering_distance_cols="euclidean",
clustering_method="average",
annotation_col = annTab, annotation_colors = anno_colors, fontsize_row=fontsize_row)
# value;
}
xx<-data.frame(group=as.factor(annTab[,cName]))
rownames(xx)<-rownames(annTab)
pheatmap(value, cluster_rows = cRow, cluster_cols = cCol,
show_colnames = FALSE, annotation_names_row = FALSE,
## annotation_col = annTab,
annotation_col = xx,
annotation_names_col = FALSE, annotation_colors = anno_colors, fontsize_row=fontsize_row)
}
#' plot gpa res
#'
#' plot gpa res
#'
#' @param gpaRes gpRes
#' @param legend whether to display it
#'
#' @return ggplot
#'
#' @export
#'
plotGPApca<-function(xres, legend=FALSE)
{
xdat<-xres$gpRes$pcaRes$pcaRes$x
grps<-xres$bundleRes$result
xdat<-xdat[names(grps),]
aDat<-data.frame(pc1=xdat[,1], pc2=xdat[,2], group=as.character(grps) )
ColorRamp <- colorRampPalette(rev(brewer.pal(n = 12,name = "Paired")))(length(unique(aDat$group)))
if(legend){
ans<-ggplot(aDat, aes(x=pc1, y=pc2, colour=group) ) + geom_point(pch=19, alpha=3/4, size=.5) + theme_bw() + scale_colour_manual(values=ColorRamp)
}
else{
ans<-ggplot(aDat, aes(x=pc1, y=pc2, colour=group) ) + geom_point(pch=19, alpha=3/4, size=.5) + theme_bw() + scale_colour_manual(values=ColorRamp) + theme(legend.position="none")
}
ans
}
#' plot tsne results
#'
#' plot tsne results
#'
#' @param mcRes, which must have TSNE.1 and TSNE.2, and group columns
#' @param legend whether to display it
#'
#' @return ggplot
#'
#' @export
#'
plotDBscan<-function(mcRes, legend=FALSE)
{
ColorRamp <- colorRampPalette(rev(brewer.pal(n = 12,name = "Paired")))(length(unique(mcRes$group)))
if(legend){
ans<-ggplot(mcRes, aes(x=TSNE.1, y=TSNE.2, colour=group) ) + geom_point(pch=19, alpha=3/4, size=.5) + theme_bw() + scale_colour_manual(values=ColorRamp)
}
else{
ans<-ggplot(mcRes, aes(x=TSNE.1, y=TSNE.2, colour=group) ) + geom_point(pch=19, alpha=3/4, size=.5) + theme_bw() + scale_colour_manual(values=ColorRamp) + theme(legend.position="none")
}
ans
}
#' @export
sc_plot_statTab<-function# multi plot of mu, alpha, mean, cov, fano, and mean vs cov
(statTab){
plot1<-ggplot(statTab, aes(x=mu)) + geom_histogram(aes(y=..density..), colour="black", fill="#238b45") + theme_bw() + ggtitle("mu")
plot2<-ggplot(statTab, aes(x=alpha)) + geom_histogram(aes(y=..density..), colour="black", fill="#88419d") + theme_bw() + ggtitle("alpha")
plot3<-ggplot(statTab, aes(x=overall_mean)) + geom_histogram(aes(y=..density..), colour="black", fill="#2b8cbe") + theme_bw() + ggtitle("overall_mean")
plot4<-ggplot(statTab, aes(x=fano)) + geom_histogram(aes(y=..density..), colour="black", fill="#d7301f") + theme_bw() + ggtitle("fano")
plot5<-ggplot(statTab, aes(x=cov)) + geom_histogram(aes(y=..density..), colour="black", fill="#feb24c") + theme_bw() + ggtitle("cov")
plot6<-ggplot(statTab, aes(x=overall_mean, y=sd)) + geom_point(pch=19, alpha=2/4, size=1, colour="#9e9ac8") + theme_bw() + ggtitle("overall mean vs sd")
multiplot(plot1, plot2, plot3, plot4, plot5, plot6, cols=2)
}
simple_tnse<-function
(tsneRes){
xres<-as.data.frame(tsneRes)
ggplot(xres, aes(x=TSNE.1, y=TSNE.2) ) + geom_point(pch=19, alpha=2/4, size=1) + theme_bw()
}
#' @export
ptsne<-function(xres, cname="study_id")
{
xi<-which(colnames(xres)==cname)
colnames(xres)[xi]<-"group"
ColorRamp <- colorRampPalette(rev(brewer.pal(n = 12,name = "Paired")))(length(unique(xres$group)))
ggplot(xres, aes(x=TSNE.1, y=TSNE.2, colour=group) ) + geom_point(pch=19, alpha=3/4, size=1) + theme_bw() + scale_colour_manual(values=ColorRamp) #+ facet_wrap( ~ k, nrow=3)
}
#' @export
plot_tsne<-function(sampTab, tsRes, cname="study_id", themeWhite=TRUE)
{
xres<-cbind(sampTab, tsRes)
xi<-which(colnames(xres)==cname)
colnames(xres)[xi]<-"group"
ColorRamp <- colorRampPalette(rev(brewer.pal(n = 12,name = "Paired")))(length(unique(xres$group)))
res<-ggplot(xres, aes(x=TSNE.1, y=TSNE.2, colour=group, fill=group) ) + geom_point(pch=21, alpha=3/4, size=1, stroke=0) + scale_colour_manual(values=ColorRamp) + scale_fill_manual(values=ColorRamp) #+ facet_wrap( ~ k, nrow=3)
if(themeWhite){
res <- res + theme_bw()
}
else{
res <- res + theme_black()
}
res
}
#' @export
tsneMult<-function# facet tsne plot by gene
(tsneDat, # cols TSNE.1, TNSE.2, and genes
genesToPlot, # genes to plot
colorPal="BuPu",
revCol=TRUE
){
require(tidyr)
tsne<-tsneDat[,c("TSNE.1", "TSNE.2", genesToPlot)]
tsneLong<-gather_(tsne, "gene", "expression", genesToPlot)
if(revCol){
ColorRamp <- rev(colorRampPalette(rev(brewer.pal(n = 7,name = colorPal)))(100))[10:100]
}
else{
ColorRamp <- colorRampPalette(rev(brewer.pal(n = 7,name = colorPal)))(100)
}
ggplot(tsneLong, aes(x=TSNE.1, y=TSNE.2, colour=expression) ) +
geom_point(pch=19, alpha=2/4, size=.25) +
theme_bw() +
scale_colour_gradientn(colours=ColorRamp) +
facet_wrap( ~ gene)
}
#' plot tsne and genes
#'
#' plot tsne and genes
#'
#' @param tsneDat tsne matrix
#' @param expDat expression matrix
#' @param genesToPlot genes to plot
#' @param colorPal "BuPu"
#' @param revCol TRUE
#' @param toScale scale?
#' @param limits
#' @return ggplot
#'
#' @export
#'
tsneMultsimp<-function(
tsneDat, # cols TSNE.1, TNSE.2
expDat,
genesToPlot, # genes to plot
colorPal="BuPu",
revCol=TRUE,
toScale=TRUE,
limits=c(-3,3),
themeWhite=TRUE
){
require(tidyr)
allgenes<-rownames(expDat)
missing<-setdiff(genesToPlot, allgenes)
cat(paste0("missing ",missing,collapse=","),"\n")
genesToPlot<-intersect(genesToPlot, allgenes)
value<-expDat[genesToPlot,]
if(toScale){
value <- t(scale(t(value)))
}
value[value < limits[1]] <- limits[1]
value[value > limits[2]] <- limits[2]
tsne<-as.data.frame(tsneDat[,c("TSNE.1", "TSNE.2")])
tsne<-cbind(tsne, t(value))
tsneLong<-gather_(tsne, "gene", "expression", genesToPlot)
tsneLong$gene_f = factor(tsneLong$gene, levels=genesToPlot)
if(revCol){
ColorRamp <- rev(colorRampPalette(rev(brewer.pal(n = 7,name = colorPal)))(100))[10:100]
}
else{
ColorRamp <- colorRampPalette(rev(brewer.pal(n = 7,name = colorPal)))(100)
}
res<- ggplot(tsneLong, aes(x=TSNE.1, y=TSNE.2, colour=expression) ) +
geom_point(pch=19, alpha=2/4, size=.25) +
scale_colour_gradientn(colours=ColorRamp) +
facet_wrap( ~ gene_f)
if(themeWhite){
res<-res + theme_bw()
}
else{
res<-res + theme_black()
}
res
# tsneLong
}
# generic red/blue heatmap
#' @export
mp_hmVars<-function# basic heatmap
(expDat,
### numeric matrix
genes,
### rownames to include in the heatmap
main='',
### optional title for teh top of the HM
clusterR=T,
### whether to cluster the Rows
clusterC=F,
### whether to cluster the columns
scale='row',
### normalize the data: 'row', 'column', or 'none'
big=FALSE,
### if not big then add cell separators
dist=utils_myDist,
###
margin=c(12,6),
###
RowSideColors=NULL,
###
ColSideColors=NULL,
###
cexCol=1,
cexRow=1,
ccol=''
){
genes<-intersect(rownames(expDat), genes);
if(length(ccol)==4){
ccol<-colorpanel(ccol$n, ccol$low,ccol$mid, ccol$high);
}
else{
ccol<-bluered(100);
}
if(is.null(RowSideColors)){
RowSideColors<-rep('white', length(genes));
}
if(is.null(ColSideColors)){
ColSideColors<-rep('white', ncol(expDat));
}
if(!clusterR & ! clusterC){
dendrogram='none';
}
if(clusterR & clusterC){
dendrogram='both';
}
if(clusterR & !clusterC){
dendrogram='row';
}
if(!clusterR & clusterC){
dendrogram='column';
}
if(!big){
heatmap.2(expDat[genes,],
#col=bluered(100),
col=ccol,
scale=scale,
trace='none',
key=T,
dendrogram=dendrogram,
Rowv=clusterR,
Colv=clusterC,
density.info='none',
margin=margin,
colsep=seq(ncol(expDat)),
rowsep=seq(length(genes)),
sepcol='white',
sepwidth=c(0.001,0.00),
main=main,
dist=dist,
RowSideColors=RowSideColors,
ColSideColors=ColSideColors,
cexCol=cexCol,
cexRow=cexRow);
}
else{
heatmap.2(expDat[genes,],col=ccol, scale=scale, trace='none', key=T,dendrogram=dendrogram,Rowv=clusterR,Colv=clusterC,density.info='none',margin=margin,main=main,dist=dist,labRow='',labCol='',RowSideColors=RowSideColors,ColSideColors=ColSideColors,cexCol=cexCol,cexRow=cexRow);
}
}
# Multiple plot function
#
# ggplot objects can be passed in ..., or to plotlist (as a list of ggplot objects)
# - cols: Number of columns in layout
# - layout: A matrix specifying the layout. If present, 'cols' is ignored.
#
# If the layout is something like matrix(c(1,2,3,3), nrow=2, byrow=TRUE),
# then plot 1 will go in the upper left, 2 will go in the upper right, and
# 3 will go all the way across the bottom.
#
# See: http://www.cookbook-r.com/Graphs/Multiple_graphs_on_one_page_(ggplot2)/
#' @export
multiplot <- function(..., plotlist=NULL, file, cols=1, layout=NULL) {
library(grid)
# Make a list from the ... arguments and plotlist
plots <- c(list(...), plotlist)
numPlots = length(plots)
# If layout is NULL, then use 'cols' to determine layout
if (is.null(layout)) {
# Make the panel
# ncol: Number of columns of plots
# nrow: Number of rows needed, calculated from # of cols
layout <- matrix(seq(1, cols * ceiling(numPlots/cols)),
ncol = cols, nrow = ceiling(numPlots/cols))
}
if (numPlots==1) {
print(plots[[1]])
} else {
# Set up the page
grid.newpage()
pushViewport(viewport(layout = grid.layout(nrow(layout), ncol(layout))))
# Make each plot, in the correct location
for (i in 1:numPlots) {
# Get the i,j matrix positions of the regions that contain this subplot
matchidx <- as.data.frame(which(layout == i, arr.ind = TRUE))
print(plots[[i]], vp = viewport(layout.pos.row = matchidx$row,
layout.pos.col = matchidx$col))
}
}
}
pcaPlot_recRes<-function(recRes, nodeName="L1_G1",legend=TRUE)
{
gpaRes<-recRes$results[[nodeName]]$gpRes
bundleRes<-recRes$results[[nodeName]]$bundleRes
aDat<-data.frame(pc1=gpaRes$pcaRes$pcaRes$x[,1], pc2=gpaRes$pcaRes$pcaRes$x[,2], group=as.character(bundleRes$result) )
ColorRamp <- colorRampPalette(rev(brewer.pal(n = 12,name = "Paired")))(length(unique(aDat$group)))
if(legend){
ans<-ggplot(aDat, aes(x=pc1, y=pc2, colour=group, pch=group) ) + geom_point(alpha=3/4, size=.5) + theme_bw() + scale_colour_manual(values=ColorRamp)
}
else{
ans<-ggplot(aDat, aes(x=pc1, y=pc2, colour=group) ) + geom_point(pch=19, alpha=3/4, size=.5) + theme_bw() + scale_colour_manual(values=ColorRamp) + theme(legend.position="none")
}
ans
}
plot_bundle<-function(gpaRes, bundleRes, legend=TRUE)
{
aDat<-data.frame(pc1=gpaRes$pcaRes$pcaRes$x[,1], pc2=gpaRes$pcaRes$pcaRes$x[,2], group=as.character(bundleRes$result) )
ColorRamp <- colorRampPalette(rev(brewer.pal(n = 12,name = "Paired")))(length(unique(aDat$group)))
if(legend){
ans<-ggplot(aDat, aes(x=pc1, y=pc2, colour=group) ) + geom_point(pch=19, alpha=3/4, size=.5) + theme_bw() + scale_colour_manual(values=ColorRamp)
}
else{
ans<-ggplot(aDat, aes(x=pc1, y=pc2, colour=group) ) + geom_point(pch=19, alpha=3/4, size=.5) + theme_bw() + scale_colour_manual(values=ColorRamp) + theme(legend.position="none")
}
ans
}
# theme_dviz_hgrid from https://github.com/clauswilke/dataviz/blob/master/R/themes.R
#' @export
theme_dviz_hgrid <- function(font_size = 14, font_family = "") {
color = "grey90"
line_size = 0.5
# Starts with theme_cowplot and then modify some parts
theme_cowplot(font_size = font_size, font_family = font_family) %+replace%
theme(
# make horizontal grid lines
panel.grid.major = element_line(colour = color,
size = line_size),
panel.grid.major.x = element_blank(),
# adjust axis tickmarks
axis.ticks = element_line(colour = color, size = line_size),
# adjust x axis
axis.line.x = element_line(colour = color, size = line_size),
# no y axis line
axis.line.y = element_blank() )
}
### theme_dviz_hgrid FROM https://github.com/clauswilke/dataviz/blob/master/R/themes.R
theme_dviz_hgrid <- function(font_size = 14, font_family = "") {
color = "grey90"
line_size = 0.5
# Starts with theme_cowplot and then modify some parts
theme_cowplot(font_size = font_size, font_family = font_family) %+replace%
theme(
# make horizontal grid lines
panel.grid.major = element_line(colour = color,
size = line_size),
panel.grid.major.x = element_blank(),
# adjust axis tickmarks
axis.ticks = element_line(colour = color, size = line_size),
# adjust x axis
axis.line.x = element_line(colour = color, size = line_size),
# no y axis line
axis.line.y = element_blank()
)
}
# theme_black from https://gist.github.com/jslefche/eff85ef06b4705e6efbc
#' @export
theme_black = function(base_size = 12, base_family = "") {
library(gridExtra)
theme_grey(base_size = base_size, base_family = base_family) %+replace%
theme(
# Specify axis options
axis.line = element_blank(),
axis.text.x = element_text(size = base_size*0.8, color = "white", lineheight = 0.9),
axis.text.y = element_text(size = base_size*0.8, color = "white", lineheight = 0.9),
axis.ticks = element_line(color = "white", size = 0.2),
#axis.title.x = element_text(size = base_size, color = "white", margin = margin(0, 10, 0, 0)),
#axis.title.y = element_text(size = base_size, color = "white", angle = 90, margin = margin(0, 10, 0, 0)),
axis.ticks.length = unit(0.3, "lines"),
# Specify legend options
legend.background = element_rect(color = NA, fill = "black"),
legend.key = element_rect(color = "white", fill = "black"),
legend.key.size = unit(1.2, "lines"),
legend.key.height = NULL,
legend.key.width = NULL,
legend.text = element_text(size = base_size*0.8, color = "white"),
legend.title = element_text(size = base_size*0.8, face = "bold", hjust = 0, color = "white"),
legend.position = "right",
legend.text.align = NULL,
legend.title.align = NULL,
legend.direction = "vertical",
legend.box = NULL,
# Specify panel options
panel.background = element_rect(fill = "black", color = NA),
panel.border = element_rect(fill = NA, color = "white"),
panel.grid.major = element_line(color = "grey35"),
panel.grid.minor = element_line(color = "grey20"),
#panel.margin = unit(0.5, "lines"),
# Specify facetting options
###strip.background = element_rect(fill = "grey30", color = "grey10"),
strip.background = element_blank(),
###strip.text.x = element_text(size = base_size*0.8, color = "white"),
strip.text.x = element_text(size = base_size, color = "white"),
strip.text.y = element_text(size = base_size*0.8, color = "white",angle = -90),
# Specify plot options
plot.background = element_rect(color = "black", fill = "black"),
###plot.title = element_text(size = base_size*1.2, color = "white"),
plot.title = element_text(size = base_size, color = "white"),
##plot.margin = unit(rep(1, 4), "lines")
)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.