|
a |
|
b/R/amaretto_vizualize.R |
|
|
1 |
#' AMARETTO_VisualizeModule |
|
|
2 |
#' |
|
|
3 |
#' Function to visualize the gene modules |
|
|
4 |
#' |
|
|
5 |
#' @param AMARETTOinit List output from AMARETTO_Initialize(). |
|
|
6 |
#' @param AMARETTOresults List output from AMARETTO_Run(). |
|
|
7 |
#' @param ProcessedData List of processed input data |
|
|
8 |
#' @param ModuleNr Module number to visualize |
|
|
9 |
#' @param SAMPLE_annotation Matrix or Dataframe with sample annotation |
|
|
10 |
#' @param ID Column used as sample name |
|
|
11 |
#' @param show_row_names If TRUE, row names will be shown on the plot. |
|
|
12 |
#' @param order_samples Order samples in heatmap by mean or by clustering |
|
|
13 |
#' |
|
|
14 |
#' @importFrom circlize colorRamp2 rand_color |
|
|
15 |
#' @importFrom grid gpar unit |
|
|
16 |
#' @importFrom stats dist hclust |
|
|
17 |
#' @importFrom dplyr left_join mutate select summarise rename filter case_when |
|
|
18 |
#' @import grDevices |
|
|
19 |
#' @import methods |
|
|
20 |
#' @importFrom ComplexHeatmap HeatmapAnnotation Heatmap draw |
|
|
21 |
#' @importFrom tibble column_to_rownames rownames_to_column |
|
|
22 |
#' @return result |
|
|
23 |
#' @export |
|
|
24 |
#' |
|
|
25 |
#' @examples |
|
|
26 |
#' data('ProcessedDataLIHC') |
|
|
27 |
#' AMARETTOinit <- AMARETTO_Initialize(ProcessedData = ProcessedDataLIHC, |
|
|
28 |
#' NrModules = 2, VarPercentage = 50) |
|
|
29 |
#' |
|
|
30 |
#' AMARETTOresults <- AMARETTO_Run(AMARETTOinit) |
|
|
31 |
#' |
|
|
32 |
#' AMARETTO_VisualizeModule(AMARETTOinit = AMARETTOinit,AMARETTOresults = AMARETTOresults, |
|
|
33 |
#' ProcessedData = ProcessedDataLIHC, ModuleNr = 1) |
|
|
34 |
AMARETTO_VisualizeModule <- function(AMARETTOinit,AMARETTOresults,ProcessedData,ModuleNr,show_row_names=FALSE, SAMPLE_annotation=NULL,ID=NULL,order_samples=NULL) { |
|
|
35 |
CNV_matrix <- ProcessedData[[2]] |
|
|
36 |
MET_matrix <- ProcessedData[[3]] |
|
|
37 |
CNVMet_Alterations <- DriversList_Alterations <- MET <- HGNC_symbol <- CNV <- NULL |
|
|
38 |
if (ModuleNr>AMARETTOresults$NrModules){ |
|
|
39 |
stop('\tCannot plot Module',ModuleNr,' since the total number of modules is',AMARETTOresults$N,'.\n') |
|
|
40 |
} |
|
|
41 |
ModuleData<-as.data.frame(AMARETTOinit$MA_matrix_Var)[AMARETTOresults$ModuleMembership==ModuleNr,] |
|
|
42 |
ModuleRegulators <- AMARETTOresults$AllRegulators[which(AMARETTOresults$RegulatoryPrograms[ModuleNr,] != 0)] |
|
|
43 |
RegulatorData <- as.data.frame(AMARETTOinit$RegulatorData)[ModuleRegulators,] |
|
|
44 |
ModuleGenes <- rownames(ModuleData) |
|
|
45 |
cat('Module',ModuleNr,'has',length(rownames(ModuleData)),'genes and',length(ModuleRegulators),'regulators for',length(colnames(ModuleData)),'samples.\n') |
|
|
46 |
Alterations<- tibble::rownames_to_column(as.data.frame(AMARETTOinit$RegulatorAlterations$Summary),"HGNC_symbol") %>% dplyr::rename(DriverList="Driver List") %>% dplyr::filter(HGNC_symbol %in% ModuleRegulators) |
|
|
47 |
Alterations<- Alterations %>% dplyr::mutate(CNVMet_Alterations=case_when(MET==1 & CNV==1~"Methylation and copy number alterations", |
|
|
48 |
CNV==1~"Copy number alterations", |
|
|
49 |
MET==1~"Methylation aberrations", |
|
|
50 |
MET==0 & CNV==0 ~"Not Altered"), |
|
|
51 |
DriversList_Alterations=case_when(DriverList==0~"Driver not predefined", |
|
|
52 |
DriverList==1~"Driver predefined")) |
|
|
53 |
|
|
|
54 |
ha_drivers <- ComplexHeatmap::HeatmapAnnotation(df = tibble::column_to_rownames(Alterations,"HGNC_symbol") %>% dplyr::select(CNVMet_Alterations,DriversList_Alterations), col = list(CNVMet_Alterations= c("Copy number alterations"="#eca400","Methylation aberrations"="#006992","Methylation and copy number alterations"="#d95d39","Not Altered"="lightgray"),DriversList_Alterations=c("Driver not predefined"="lightgray","Driver predefined"="#588B5B")),which = "column", height = grid::unit(0.3, "cm"),name="", |
|
|
55 |
annotation_legend_param = list(title_gp = grid::gpar(fontsize = 8),labels_gp = grid::gpar(fontsize = 6))) |
|
|
56 |
|
|
|
57 |
overlapping_samples <- colnames(ModuleData) |
|
|
58 |
|
|
|
59 |
# Add NAs for Gene Expression samples not existing in CNV or MET data |
|
|
60 |
|
|
|
61 |
if (!is.null(CNV_matrix)){ |
|
|
62 |
non_CNV_sample_names<- overlapping_samples[!overlapping_samples%in%colnames(CNV_matrix)] |
|
|
63 |
print(non_CNV_sample_names) |
|
|
64 |
print(nrow(CNV_matrix)) |
|
|
65 |
print(length(non_CNV_sample_names)) |
|
|
66 |
CNV_extension_mat<- matrix(data=NA,nrow=nrow(CNV_matrix),ncol=length(non_CNV_sample_names)) |
|
|
67 |
colnames(CNV_extension_mat)<-non_CNV_sample_names |
|
|
68 |
CNV_matrix<-cbind(CNV_matrix,CNV_extension_mat) |
|
|
69 |
} |
|
|
70 |
if (!is.null(MET_matrix)){ |
|
|
71 |
non_MET_sample_names<- overlapping_samples[!overlapping_samples%in%colnames(MET_matrix)] |
|
|
72 |
MET_extension_mat<- matrix(data=NA,nrow=nrow(MET_matrix),ncol=length(non_MET_sample_names)) |
|
|
73 |
colnames(MET_extension_mat)<-non_MET_sample_names |
|
|
74 |
MET_matrix<-cbind(MET_matrix,MET_extension_mat) |
|
|
75 |
} |
|
|
76 |
|
|
|
77 |
if(is.null(order_samples)){ |
|
|
78 |
overlapping_samples_clust<-overlapping_samples[order(colMeans(ModuleData[,overlapping_samples]))] |
|
|
79 |
}else if(order_samples=="clust"){ |
|
|
80 |
SampleClustering<-stats::hclust(stats::dist(t(ModuleData[,overlapping_samples])), method = "complete", members = NULL) |
|
|
81 |
overlapping_samples_clust<-overlapping_samples[SampleClustering$order] |
|
|
82 |
}else { |
|
|
83 |
print("ordering type not recognized, samples will be orderd based on mean expression of the module genes") |
|
|
84 |
overlapping_samples_clust<-overlapping_samples[order(colMeans(ModuleData[,overlapping_samples]))] |
|
|
85 |
} |
|
|
86 |
|
|
|
87 |
ClustRegulatorData <- t(RegulatorData[,overlapping_samples_clust]) |
|
|
88 |
ClustModuleData <- t(ModuleData[,overlapping_samples_clust]) |
|
|
89 |
|
|
|
90 |
if(length(ClustModuleData)<50){ |
|
|
91 |
fontsizeMo=6 |
|
|
92 |
} else if (length(ClustModuleData)<200){ |
|
|
93 |
fontsizeMo=4 |
|
|
94 |
} else {fontsizeMo=2} |
|
|
95 |
|
|
|
96 |
Regwidth <- ncol(ClustRegulatorData)*0.5 |
|
|
97 |
ha_Reg <- Heatmap(ClustRegulatorData, name = "Gene Expression", column_title = "Regulator Genes\nExpression",cluster_rows=FALSE,cluster_columns=TRUE,show_column_dend=FALSE,show_column_names=TRUE,show_row_names=FALSE,column_names_gp = grid::gpar(fontsize = fontsizeMo),top_annotation = ha_drivers, |
|
|
98 |
column_title_gp = grid::gpar(fontsize = 6, fontface = "bold"), col=circlize::colorRamp2(c(-max(abs(ClustRegulatorData)), 0, max(abs(ClustRegulatorData))), c("darkblue", "white", "darkred")),heatmap_legend_param = list(color_bar = "continuous",legend_direction = "horizontal",title_gp = grid::gpar(fontsize = 8),labels_gp = grid::gpar(fontsize = 6)), width = grid::unit(Regwidth, "cm")) |
|
|
99 |
|
|
|
100 |
|
|
|
101 |
ha_Module <- Heatmap(ClustModuleData, name = " ", column_title = "Target Genes\nExpression",cluster_rows=FALSE,cluster_columns=TRUE,show_column_dend=FALSE,show_column_names=TRUE,show_row_names=show_row_names,column_names_gp = grid::gpar(fontsize = fontsizeMo),show_heatmap_legend = FALSE, |
|
|
102 |
column_title_gp = grid::gpar(fontsize = 6, fontface = "bold"), col=circlize::colorRamp2(c(-max(abs(ClustModuleData)), 0, max(abs(ClustModuleData))), c("darkblue", "white", "darkred")),heatmap_legend_param = list(color_bar = "continuous",legend_direction = "horizontal",title_gp = grid::gpar(fontsize = 8),labels_gp = grid::gpar(fontsize = 6))) |
|
|
103 |
|
|
|
104 |
ha_list<- ha_Reg + ha_Module |
|
|
105 |
if (!is.null(MET_matrix)){ |
|
|
106 |
METreg <- intersect(rownames(AMARETTOinit$RegulatorAlterations$MET),ModuleRegulators) |
|
|
107 |
print("MET regulators will be included when available") |
|
|
108 |
if (length(METreg)>0){ |
|
|
109 |
MET_matrix = as.data.frame(MET_matrix) |
|
|
110 |
METData2 = METData = as.matrix(MET_matrix[unlist(Alterations %>% dplyr::filter(MET==1) %>% dplyr::select(HGNC_symbol)),overlapping_samples_clust]) |
|
|
111 |
METData2[which(METData>0)] <- "Hyper-methylated" # hyper |
|
|
112 |
METData2[which(METData<0)] <- "Hypo-methylated" # hypo |
|
|
113 |
METData2[which(METData==0)] <- "Not altered" # nothing |
|
|
114 |
METData2<-t(METData2) |
|
|
115 |
Metwidth=ncol(METData2)*0.7 |
|
|
116 |
Met_col=structure(c("#006992","#d95d39","white"),names=c("Hyper-methylated","Hypo-methylated","Not altered")) |
|
|
117 |
ha_Met <- Heatmap(METData2, name = "Methylation State", column_title = "Methylation", cluster_rows=FALSE,cluster_columns=TRUE,show_column_dend=FALSE,show_column_names=TRUE,show_row_names=FALSE,column_names_gp = grid::gpar(fontsize = 6),show_heatmap_legend = TRUE, |
|
|
118 |
column_title_gp = gpar(fontsize = 6, fontface = "bold"), col = Met_col, width = grid::unit(Metwidth, "cm"),heatmap_legend_param = list(title_gp = gpar(fontsize = 8),labels_gp = grid::gpar(fontsize = 6))) |
|
|
119 |
ha_list<- ha_Met + ha_list |
|
|
120 |
} |
|
|
121 |
} |
|
|
122 |
|
|
|
123 |
if (!is.null(CNV_matrix)){ |
|
|
124 |
CNVreg <- intersect(rownames(AMARETTOinit$RegulatorAlterations$CNV),ModuleRegulators) |
|
|
125 |
print("CNV regulators will be included when available") |
|
|
126 |
if (length(CNVreg)>0){ |
|
|
127 |
CNV_matrix = as.data.frame(CNV_matrix) |
|
|
128 |
CNVData2 = CNVData = as.matrix(CNV_matrix[unlist(Alterations %>% dplyr::filter(CNV==1) %>% dplyr::select(HGNC_symbol)),overlapping_samples_clust]) |
|
|
129 |
CNVData2[which(CNVData>=0.1)] <- "Amplified" # amplification |
|
|
130 |
CNVData2[which(CNVData<=(-0.1))] <- "Deleted" # deletion |
|
|
131 |
CNVData2[which(CNVData<0.1 & CNVData>(-0.1))] <- "Not altered" # nothing |
|
|
132 |
CNVData2<-t(CNVData2) |
|
|
133 |
CNVwidth=ncol(CNVData2)*0.7 |
|
|
134 |
CNV_col=structure(c("#006992","#d95d39","white"),names=c("Deleted","Amplified","Not altered")) |
|
|
135 |
ha_CNV <- Heatmap(CNVData2, name = "CNV State", column_title = "CNV", cluster_rows=FALSE,cluster_columns=TRUE,show_column_dend=FALSE,show_column_names=TRUE,show_row_names=FALSE,column_names_gp = grid::gpar(fontsize = 6),show_heatmap_legend = TRUE, |
|
|
136 |
column_title_gp = grid::gpar(fontsize = 6, fontface = "bold"),col = CNV_col,width = grid::unit(CNVwidth, "cm"),heatmap_legend_param = list(title_gp = grid::gpar(fontsize = 8),labels_gp = grid::gpar(fontsize = 6))) |
|
|
137 |
ha_list<-ha_CNV + ha_list |
|
|
138 |
} |
|
|
139 |
} |
|
|
140 |
|
|
|
141 |
if (!is.null(SAMPLE_annotation)){ |
|
|
142 |
if (ID %in% colnames(SAMPLE_annotation)){ |
|
|
143 |
SAMPLE_annotation_fil<-as.data.frame(SAMPLE_annotation) %>% dplyr::filter(!!as.name(ID) %in% overlapping_samples_clust) |
|
|
144 |
suppressWarnings(SAMPLE_annotation_fil<-dplyr::left_join(as.data.frame(overlapping_samples_clust),SAMPLE_annotation_fil,by=c("overlapping_samples_clust"=ID))) |
|
|
145 |
SAMPLE_annotation_fil<-tibble::column_to_rownames(SAMPLE_annotation_fil,"overlapping_samples_clust") |
|
|
146 |
cat(nrow(SAMPLE_annotation_fil),"samples have an annotation.\n") |
|
|
147 |
cat(ncol(SAMPLE_annotation_fil),"annotations are added") |
|
|
148 |
#define colors |
|
|
149 |
col<-c() |
|
|
150 |
for (sample_column in colnames(SAMPLE_annotation_fil)[colnames(SAMPLE_annotation_fil) != ID]){ |
|
|
151 |
newcol<-circlize::rand_color(n=length(unique(SAMPLE_annotation_fil[,sample_column])),luminosity = "bright") |
|
|
152 |
names(newcol)<-unique(SAMPLE_annotation_fil[,sample_column]) |
|
|
153 |
col<-c(col,newcol) |
|
|
154 |
} |
|
|
155 |
ha_anot<-Heatmap(SAMPLE_annotation_fil, name="Sample Annotation", column_title = "Sample\nAnnotation", column_title_gp = grid::gpar(fontsize = 6, fontface = "bold"), col=col, show_row_names=FALSE,width = unit(ncol(SAMPLE_annotation_fil) * 2, "mm"), |
|
|
156 |
column_names_gp = gpar(fontsize = 6),heatmap_legend_param = list(title_gp = grid::gpar(fontsize = 8),labels_gp = grid::gpar(fontsize = 6))) |
|
|
157 |
ha_list<-ha_list + ha_anot |
|
|
158 |
} else {print("The ID is not identified as a column name in the annotation")} |
|
|
159 |
} |
|
|
160 |
ComplexHeatmap::draw(ha_list,heatmap_legend_side = "bottom",annotation_legend_side="bottom") |
|
|
161 |
} |