gravatar for Eugene A

2 hours ago by

These are from R pheatmap package (www.rdocumentation.org/packages/pheatmap/versions/1.0.12/topics/pheatmap, davetang.org/muse/2018/05/15/making-a-heatmap-in-r-with-the-pheatmap-package/)

the code will be according to the following line:

library(pheatmap)
library(RColorBrewer)
plot_heatmap <- function(gene_counts, gene_set, gene_ann = NA, samp_ann = NA, change_rows = FALSE){
  #function to plot heatmap with dendrograms
  #input - dataframe with normalized counts (columns - samples, rows - genes) and set of genes (vector) to perform clasterisation on
  # gene_ann = NA, samp_ann = NA 1 column dataframes for legends (side-bars)
  #1st - prefrom clast on samples based on cpearman correlation
  #2nd - perform clast on genes based on spearman correlation
  #3rd - draw resulting heatmap of log(gene expression), ordering genes and samples based on previouse clasterisaiton and adding previouse dendros

  breaksList <- seq(-1, 1, by = 0.1)

  #class names for legends and side-bars
  if(!is.na(gene_ann)[1]){colnames(gene_ann) <- c('Genes')}
  if(!is.na(samp_ann)[1]){colnames(samp_ann) <- c('Sample_types')}

  # get necessary maatrix
  gene_mtr <- gene_counts[gene_set,]

  if(change_rows){
    colnames(gene_mtr) <- row.names(samp_ann)
  }
  #clustering of samples based on correlation
  correlationmatrix <- cor(gene_mtr, use="pairwise.complete.obs", method = 'spearman')
  samp_pm <- pheatmap(correlationmatrix, main = 'Spearman', silent = T)
  samp_order <- samp_pm$tree_row$order
  #clustering of genes based on correlation
  correlationmatrix <- cor(t(gene_mtr), use="pairwise.complete.obs", method = 'spearman')
  gene_pm <- pheatmap(correlationmatrix, main = 'Spearman', silent = T)
  gene_order <- gene_pm$tree_row$order

  #make a heatmap of log(expression) with previously identified clusters
  hetM <- pheatmap(log(gene_mtr+1), cluster_rows = gene_pm$tree_col, cluster_cols = samp_pm$tree_row,
                   #graffical parameters go here
                   color = colorRampPalette(rev(brewer.pal(n = 7, name = "RdYlBu")))(length(breaksList)),
                   annotation_col = samp_ann, fontsize =7,
                   annotation_row = gene_ann
  )

  #return(hetM)
}

Best



Source link