Commit c78a7f29 authored by smorabit's avatar smorabit
Browse files

kME function and plotting function

parent 1f6871b7
Loading
Loading
Loading
Loading
+50 −1
Original line number Diff line number Diff line
@@ -380,6 +380,24 @@ ModuleEigengenes <- function(seurat_obj, group.by.vars=NULL, verbose=TRUE, ...){
  seurat_obj
}

#' GetMEs
#'
#' Function to retrieve module eigengens from Seurat object.
#'
#' @param seurat_obj A Seurat object
#' @keywords scRNA-seq
#' @export
#' @examples
#'  ModuleEigengenes(pbmc)
GetMEs <- function(seurat_obj, harmonized=TRUE){
  if(harmonized == TRUE && !is.null(seurat_obj@misc$wgcna_hMEs)){
    MEs <- seurat_obj@misc$wgcna_hMEs
  } else{
    MEs <- seurat_obj@misc$wgcna_MEs
  }
  MEs
}


#' ModuleConnectivity
#'
@@ -390,6 +408,37 @@ ModuleEigengenes <- function(seurat_obj, group.by.vars=NULL, verbose=TRUE, ...){
#' @export
#' @examples
#' ModuleConnectivity(pbmc)
ModuleConnectivity <- function(seurat_obj, group.by.vars=NULL, verbose=TRUE, ...){
ModuleConnectivity <- function(seurat_obj, harmonized=TRUE, ...){

  # get expression matrix
  genes_use <- names(seurat_obj@misc$wgcna_net$colors)

  # datExpr for full expression dataset
  # transpose expression data and subset by genes used for WGCNA
  datExpr <- t(GetAssayData(
    seurat_obj,
    assay=seurat_obj@misc$wgcna_params$metacell_assay,
    slot=seurat_obj@misc$wgcna_params$metacell_slot
  ))[,genes_use]

  # get MEs:
  MEs <- GetMEs(seurat_obj=seurat_obj, harmonized=harmonized)

  tic("SignedKME")
  kMEs <- signedKME(
    datExpr,
    MEs,
    outputColumnName = "kME",
    corFnc = "bicor",
    ...
  )
  toc()

  # add module color to the kMEs table
  kMEs <- cbind(cur_seurat@misc$wgcna_net$colors, kMEs)
  colnames(kMEs) <- c('module', colnames(MEs))

  seurat_obj@misc$wgcna_kMEs <- kMEs
  seurat_obj

}

R/plotting.R

0 → 100644
+109 −0
Original line number Diff line number Diff line

#' PlotDendrogram
#'
#' Plot WGCNA dendrogram
#'
#' @param seurat_obj A Seurat object
#' @keywords scRNA-seq
#' @export
#' @examples
#' PlotDendrogram
PlotDendrogram <- function(
  seurat_obj, groupLabels="Module colors",
  dendroLabels = FALSE, hang = 0.03, addGuide = TRUE, guideHang = 0.05,
  main = "", ...
){

  WGCNA::plotDendroAndColors(
    seurat_obj@misc$wgcna_net$dendrograms[[1]],
    as.character(seurat_obj@misc$wgcna_net$colors),
    groupLabels=groupLabels,
    dendroLabels = dendroLabels,
    hang = hang,
    addGuide = addGuide,
    guideHang = guideHang,
    main = main,
    ...
  )
}


#' MEFeaturePlot
#'
#' Plot module eigengenes as a FeaturePlot
#'
#' @param seurat_obj A Seurat object
#' @keywords scRNA-seq
#' @export
#' @examples
#' MEFeaturePlot
MEFeaturePlot<- function(
  seurat_obj, modules,
  harmonized=TRUE, reduction='umap',
  order=TRUE, restrict_range=TRUE, point_size = 0.5,
  label_legend = FALSE
){

  # get MEs from seurat object
  MEs <- GetMEs(seurat_obj, harmonized)

  # get  reduction from seurat obj
  umap <- seurat_obj@reductions[[reduction]]@cell.embeddings
  x_name <- colnames(umap)[1]
  y_name <- colnames(umap)[2]

  # merge into one df for plotting
  plot_df <- cbind(umap, MEs) %>% as.data.frame()

  plot_list <- list()
  for(cur_mod in modules){

    cur_color <- cur_mod

    # reset the range of the plot:
    plot_range <- plot_df[,cur_mod] %>% range
    if(restrict_range){
      if(abs(plot_range[1]) > abs(plot_range[2])){
        plot_range[1] <- -1*plot_range[2]
      } else{
        plot_range[2] <- -1*plot_range[1]
      }
      plot_df[,cur_mod] <- ifelse(plot_df[,cur_mod] > plot_range[2], plot_range[2], plot_df[,cur_mod])
      plot_df[,cur_mod] <- ifelse(plot_df[,cur_mod] < plot_range[1], plot_range[1], plot_df[,cur_mod])
    }

    # order points:
    if(order){
      plot_df <- plot_df %>% dplyr::arrange_(cur_mod)
    }

    # label for plot:
    if(harmonized){
      label <- paste0('hME_', cur_mod)
    } else{
      label <- paste0('ME_', cur_mod)
    }

    # plot with ggplot
    plot_list[[cur_mod]] <- plot_df %>%
      ggplot(aes_string(x=x_name, y=y_name, color=cur_mod)) +
      geom_point(size=0.5) +
      scale_color_gradient2(
        low='grey75', mid='grey95', high=cur_color,
        breaks = c(plot_range[1], 0, plot_range[2]),
        labels = c('-', '0', '+'),
      ) +
      ggtitle(label) + umap_theme +
      labs(color="")
  }

  # return plot
  if(length(plot_list) == 1){
    p <- plot_list[[1]]
  } else{
    p <- plot_list
  }

  p

}