R/MDSPlot.R

In biomics-pasteur-fr/RNADiff: R Tools for RNA-Seq Projects at Institut Pasteur (Biomics Platform)

#' MDS plot of samples (for edgeR objects)
#'
#' MDS plot of samples (for edgeR objects)
#'
#' @param dge an \code{DGEList} object
#' @param group vector of the condition from which each sample belongs
#' @param n number of features to keep among the most variant
#' @param gene.selection \code{"pairwise"} to choose the top features separately for each pairwise
#'  comparison between the samples or \code{"common"} to select the same features for all
#'  comparisons. Only used when \code{method="logFC"}
#' @param out \code{TRUE} to export the figure
#' @param col colors of the bars
#' @param versionName versionName of the project
#' @return A MDS plot
#' @author Marie-Agnes Dillies and Hugo Varet

# created Nov 14th, 2013
# modified Mar 21st, 2014 (removed outputfile argument)
# modified Aug 5th, 2014 (removed graphDir argument)
# modified Dec 02nd, 2014 (removed dev.new() et dev.off() before and after plotMDS())
# modified June 15th, 2015 (take the mininimum between 500 and the number of features)
# modified August 28th, 2019 (ggplot2)

MDSPlot <- function(dge, group, n=min(500,nrow(dge$counts)), gene.selection=c("pairwise", "common"),
                    out = TRUE, col=c("lightblue","orange","MediumVioletRed","SpringGreen"),
                    versionName="."){
  gene.selection <- gene.selection[1]
  group <- data.frame(group=apply(group, 1, paste, collapse="-"))
  group$group <- factor(group$group, levels=unique(group$group))
  if (out) pdf(file=paste0("figures/", versionName, "-MDS.pdf"))
  coord <- plotMDS(dge, top=n, method="logFC", gene.selection=gene.selection, plot=FALSE)
  coord <- as.data.frame(coord)
  d <- data.frame(coord[,c("x", "y")], 
                  group=group, 
                  sample=factor(row.names(coord), levels=row.names(coord)))
  print(ggplot(data=d, aes(x=.data$x, y=.data$y, color=group, label=sample)) + 
          geom_point(show.legend=TRUE, size=3) +
          labs(color="") +
          scale_colour_manual(values=col) +
          geom_text_repel(show.legend=FALSE, size=5, point.padding=0.2) +
          xlab("Leading logFC dimension 1") +
          ylab("Leading logFC dimension 2") +
          ggtitle(paste(versionName, "Multi-Dimensional Scaling", sep=" - ")))
  if (out) dev.off()      
}