R/plotCorr.R

In eifer4/stochasticSampling: Stochastic partitioning sampler.

plotCorr <- function(X_corr=NULL, Y_corr = NULL, Y_grp = NULL, X_grp = NULL, file = NULL, x.title = NULL, xlab = NULL, ylab = NULL, y.title = NULL, x.breaks=NULL, y.breaks=NULL,y.text = NULL, x.text=NULL, ...) {
  
  dots <- list(...)
  
  Y.order <- X.order <- NULL
  
  if("X.order" %in% names(dots)) X.order <- dots$X.order
  if("Y.order" %in% names(dots)) Y.order <- dots$Y.order
  
  if(!is.null(X_corr)){
    
    Q <- nrow(X_corr)
    P <- ncol(X_corr)
    
    
    if(is.null(X_grp) | is.null(Y_grp)) {
      
      Xcorrdf <- data.frame(X_corr)
      
      if(is.null(Y.order)) Y.order <- rownames(Xcorrdf)
      Xcorrdf$Yid <- Y.order
      
      Xcorrdf <- melt(Xcorrdf, id.vars="Yid")
      
      if(is.null(X.order)) X.order <- as.numeric(Xcorrdf$variable)
      levels(Xcorrdf$variable) <- X.order
      
      
      Xcorrdf$Yid <- factor(Xcorrdf$Yid, levels=sort(unique(Y.order)))
      Xcorrdf$variable <- factor(Xcorrdf$variable, levels=sort(unique(X.order)))
      Xcorrdf$alpha <- 1
      breaks<- seq(0,1, 0.2)
      break.lab <- format(breaks)
      gradient <- scale_fill_gradient2(low = "#ffeda0" , mid="#feb24c", high = "#f03b20",  midpoint = 0.5, name="Posterior Probability\nof Relationship", breaks=breaks, labels=break.lab, limits=c(0,1), expand=c(0,0))
      height <- 4
      width <- 5
    } else {
      stopifnot(length(X_grp) == length(Y_grp))
      K <- length(Y_grp)
      
      grps <- matrix(1, nrow=nrow(X_corr), ncol=ncol(X_corr))
      for(i in 1:K) grps[Y_grp[[i]], X_grp[[i]]] <- 2
      
      X_corrtemp <- X_corr
      X_corrtemp[which(grps==1,arr.ind=FALSE)] <- X_corrtemp[which(grps==1,arr.ind=FALSE)]*-1
      
      grp.col <- matrix(1, nrow=nrow(X_corr), ncol=ncol(X_corr))
      grp.col[which(grps==2 & X_corrtemp<0.25 & X_corrtemp>=0,arr.ind=FALSE)] <-  0
      
      
      Xcorrdf <- data.frame(X_corrtemp)
      Xcorrdf$Yid <- rownames(Xcorrdf)
      
      Xcorrdf <- melt(Xcorrdf, id.vars="Yid")
  
      
      Xcorrdf$group <- factor(c(grps), labels=c("out", "in"))
      Xcorrdf$Yid <- factor(Xcorrdf$Yid, levels=unlist(Y_grp)[length(unlist(Y_grp)):1])
      Xcorrdf$variable <- factor(as.numeric(Xcorrdf$variable), levels=1:P)
      Xcorrdf$alpha <- c(grp.col)
      breaks<- seq(-1,1, 0.25)
      break.lab <- format(breaks)
      gradient <- scale_fill_gradient2(low = "#2b83ba" , mid="#ffffbf", high = "#d7191c",  midpoint = 0, guide='none',breaks=breaks, labels=break.lab, limits=c(-1,1), expand=c(0,0))
      height<-4.2
      width <- 4.2
    }
    
    if(is.null(x.title)) x.title <- "Correlation of Y and X Variables"
    if(is.null(ylab)) ylab <- "Y Variables"
    if(is.null(xlab)) xlab <- "X Variables"
    
    Xcorrp <- ggplot(Xcorrdf, aes(x=variable, y=Yid, fill=value, alpha=alpha)) + 
      geom_tile() + 
      theme(#axis.text=element_blank(),
            axis.ticks=element_blank(),
            panel.grid.major = element_blank(), 
            panel.grid.minor = element_blank()) +
      scale_x_discrete(breaks=x.breaks, labels=x.text, expand = c(0, 0)) +
      scale_y_discrete(breaks=y.breaks, labels=y.text, expand = c(0, 0)) +
      ggtitle(x.title) + ylab(ylab)  +
      xlab(xlab)
    Xcorrp <- Xcorrp + gradient + scale_alpha(guide="none")
    if(is.null(file)){
      print(Xcorrp)
    } else {
      pdf(file=paste0(file,"_XCorr.pdf"), height=height, width=width)
      print(Xcorrp)
      dev.off()
    }
  }
  
  if(!is.null(Y_corr)) {
    if(is.null(Y_grp)) {

      Q <- nrow(Y_corr)
      
      Ycorrdf <- data.frame(Y_corr)
      if(is.null(Y.order)){
        Y.order <- order(rownames(Ycorrdf))
      }

      Ycorrdf$Yid <- Y.order
      
      Ycorrdf <- melt(Ycorrdf, id.vars="Yid")
      levels(Ycorrdf$variable) <- Y.order

      Ycorrdf$Yid <- factor(Ycorrdf$Yid, levels=sort(unique(Y.order)))
      Ycorrdf$variable <- factor(Ycorrdf$variable, levels=sort(unique(Y.order)))
      Ycorrdf$alpha <- 1
      breaks<- seq(0,1, 0.2)
      break.lab <- format(breaks)
      gradient <- scale_fill_gradient2(low = "#ffeda0" , mid="#feb24c", high = "#f03b20",  midpoint = 0.5, name="Posterior Probability\nof Relationship", breaks=breaks, labels=break.lab, limits=c(0,1), expand=c(0,0))
    } else {
      
      K <- length(Y_grp)
      Q <- nrow(Y_corr)

      grps <- matrix(1, nrow=nrow(Y_corr), ncol=ncol(Y_corr))
      for(i in 1:10) grps[Y_grp[[i]], Y_grp[[i]]] <- 2
      
      Ycorrstemp <- Y_corr
      Ycorrstemp[which(grps==1,arr.ind=FALSE)] <- Ycorrstemp[which(grps==1,arr.ind=FALSE)]*-1
      grp.col <- matrix(1, nrow=nrow(Y_corr), ncol=ncol(Y_corr))
      grp.col[which(grps==2 & Ycorrstemp < 0.25 & Ycorrstemp>=0,arr.ind=FALSE)] <-  0
      Ycorrdf <- data.frame(Ycorrstemp)
      Ycorrdf$Yid <- rownames(Ycorrdf)
      
      Ycorrdf <- melt(Ycorrdf, id.vars="Yid")
      Yorder <- unlist(Y_grp)[length(unlist(Y_grp)):1]
      
      Ycorrdf$group <- factor(c(grps), labels=c("out", "in"))
      Ycorrdf$Yid <- factor(Ycorrdf$Yid, levels=Yorder)
      Ycorrdf$variable <- factor(as.numeric(Ycorrdf$variable), levels = unlist(Y_grp))
      Ycorrdf$alpha <- c(grp.col)
      breaks<- seq(-1,1, 0.25)
      break.lab <- format(breaks)
      gradient <- scale_fill_gradient2(low = "#2b83ba" , mid="#ffffbf", high = "#d7191c",  midpoint = 0, name="Inclusion Probability", breaks=breaks, labels=break.lab, limits=c(-1,1), expand=c(0,0))
    }
    
    if(is.null(y.title)) y.title <- "Correlation of Y Variables"
    Ycorrp <- ggplot(Ycorrdf, aes(x=variable, y=Yid, fill=value, alpha=alpha)) + 
      geom_tile() + 
      theme(#axis.text=element_blank(),
        axis.ticks=element_blank(),
        panel.grid.major = element_blank(), 
        panel.grid.minor = element_blank()) +
      scale_x_discrete(breaks=c(0,y.breaks[2]), labels=y.text, expand = c(0, 0)) +
      scale_y_discrete(breaks=y.breaks, labels=y.text, expand = c(0, 0)) +
      ylab(NULL) + xlab(NULL)+
      ggtitle(y.title)
    Ycorrp <- Ycorrp + gradient + scale_alpha(guide="none")
    if(is.null(file)){
      print(Ycorrp)
    } else {
      pdf(file=paste0(file,"_YCorr.pdf"), height=4, width=5)
      print(Ycorrp)
      dev.off()
    }
    
  }
  
}