Melanoma KMT2D H3K4me1 heatmap & average intensity curves

Mayinuer(Mahinur) Maitituoheti August 2019 based on code from Ming Tang

Load necessary programs

library(rtracklayer)
library(ComplexHeatmap)
library(EnrichedHeatmap)
library(circlize)
library(magrittr)
library(dplyr)
library(lazyeval)
library(tidyr)
library(ggplot2)

Import and normalization of ChIP-seq reads(RPKM)

# read in merged H3K4me1 bed file from KMT2D WT and Mut samples
H3K4me1_bed<- import("/Volumes/Mahinur/SKCM_KMT2D/SKCM_KMT2D_H3K4me1.merged.bed", format = "BED")
# set peaks to 10kb window centered on the middle 
H3K4me1.5kb<- resize(H3K4me1_bed, width = 10000, fix = "center")
# only read in a subset of the bigwig files using which argument!
IA99.H3K4me1<- import("/Volumes/Mahinur/SKCM_KMT2D/07bigwig/SKCM-IA99_A.bw", format = "bigWig", which = H3K4me1.5kb)
IE17.H3K4me1<- import("/Volumes/Mahinur/SKCM_KMT2D/07bigwig/SKCM-IE17_A.bw", format = "bigWig", which = H3K4me1.5kb)
# set the center for resized peaks
H3K4me1.center<- resize(H3K4me1.5kb, width =1, fix = "center")
# normalization for resized H3K4me1 ChIP-seq reads
IA99.mat1<- normalizeToMatrix(IA99.H3K4me1, H3K4me1.center, value_column = "score",
                              mean_mode="w0", w=100, extend = 5000)
IE17.mat1<- normalizeToMatrix(IE17.H3K4me1, H3K4me1.center, value_column = "score",
                              mean_mode="w0", w=100, extend = 5000)
# check quantile range outliers for the matrix 
quantile(IA99.mat1, probs = c(0.005, 0.5,0.995))
##     0.5%      50%    99.5% 
##  0.00000 11.25800 36.00915
quantile(IE17.mat1, probs = c(0.005, 0.5,0.995))
##    0.5%     50%   99.5% 
##  0.0000  9.7348 28.3125

H3K4me1 heatmap

# mapping haetmap colors and set the range
col_fun1<- circlize::colorRamp2(c(0, 40), c("white", "red"))
# set haetmap parameters
ht_global_opt(heatmap_column_names_gp = gpar(fontsize = 14),
              heatmap_legend_title_gp = gpar(fontsize = 14),
              heatmap_legend_labels_gp = gpar(fontsize = 14),
              ADD = TRUE)
# Print the heatmap
print(EnrichedHeatmap(IA99.mat1, axis_name_rot = 45, name = "KMT2D WT H3K4me1",
                column_title = "KMT2D WT", use_raster = TRUE,
                col= col_fun1, pos_line =F,
                axis_name_gp = gpar(fontsize = 14), axis_name = c("-5kb", "center", "5kb")) +
EnrichedHeatmap(IE17.mat1, axis_name_rot = 45, name = "KMT2D Mut H3K4me1",
                  column_title ="KMT2D Mut", use_raster = TRUE, col= col_fun1, pos_line =F,
                  axis_name_gp = gpar(fontsize = 14), axis_name = c("-5kb", "center", "5kb")))

Calculate average intensity of the H3K4me1 signal

# mean values of normalized H3K4me1 signal
IA99_mean<- as.data.frame(colMeans(IA99.mat1))
IE17_mean<- as.data.frame(colMeans(IE17.mat1))
# gather all samples average intensity values
H3K4me1_all<- cbind(IA99_mean, IE17_mean) %>%
  mutate(pos = rownames(IA99_mean)) %>%
  dplyr::rename(KMT2D_WT = `colMeans(IA99.mat1)`, KMT2D_Mut = `colMeans(IE17.mat1)`) %>%
  gather(sample, value, 1:2)
# set line plot position
H3K4me1_all$pos<- factor(H3K4me1_all$pos, levels= rownames(IA99_mean))

Average intensity curves for normalized H3K4me1 signal

print(ggplot(H3K4me1_all, aes(x = pos,y = value, group = sample)) + geom_line(aes(color = sample), size = 1.2) +
  theme_classic(base_size = 18) +
  theme(axis.ticks.x = element_blank()) +
  scale_x_discrete(breaks = c("u1", "d1", "d50"), labels =c ("-5kb", "center", "5kb")) +
  xlab(NULL) +
  ggtitle("KMT2D H3K4me1 signal") +
theme(panel.spacing = unit(2, "lines")))

Note that the echo = FALSE parameter was added to the code chunk to prevent printing of the R code that generated the plot.