使用ggplot2将显着性水平添加到矩阵相关热图

2 回答

• 28

  这只是为了提升最终解决方案的尝试，我在这里绘制了星星作为解决方案的指标，但正如我所说，目标是找到一个比星星更好的图形解决方案 . 我只是使用geom_point和alpha来表示显着性水平，但是NA（包括非重要值）的问题将显示为三星级重要性的问题，如何解决这个问题？我认为使用一种颜色在使用多种颜色时可能会更加眼睛友好，并避免给情节加重许多细节以便眼睛解决 . 提前致谢 .
  这是我第一次尝试的情节：

  

  或者可能这样更好？！

  

  我认为现在最好的是下面的那个，直到你想出更好的东西！
  

  根据要求，以下代码适用于最后一个热图：

  # Function to get the probability into a whole matrix not half, here is Spearman you can change it to Kendall or Pearson
  cor.prob.all <- function (X, dfr = nrow(X) - 2) {
  R <- cor(X, use="pairwise.complete.obs",method="spearman")
  r2 <- R^2
  Fstat <- r2 * dfr/(1 - r2)
  R<- 1 - pf(Fstat, 1, dfr)
  R[row(R) == col(R)] <- NA
  R
  }
  # Change matrices to dataframes
  nbar<- as.data.frame(cor(nba[2:ncol(nba)]),method="spearman") # to a dataframe for r^2
  nbap<- as.data.frame(cor.prob.all(nba[2:ncol(nba)])) # to a dataframe for p values
  # Reset rownames
  nbar <- data.frame(row=rownames(nbar),nbar) # create a column called "row" 
  rownames(nbar) <- NULL
  nbap <- data.frame(row=rownames(nbap),nbap) # create a column called "row" 
  rownames(nbap) <- NULL
  # Melt
  nbar.m <- melt(nbar)
  nbap.m <- melt(nbap)
  # Classify (you can classify differently for nbar and for nbap also)         
  nbar.m$value2<-cut(nbar.m$value,breaks=c(-1,-0.75,-0.5,-0.25,0,0.25,0.5,0.75,1),include.lowest=TRUE, label=c("(-0.75,-1)","(-0.5,-0.75)","(-0.25,-0.5)","(0,-0.25)","(0,0.25)","(0.25,0.5)","(0.5,0.75)","(0.75,1)")) # the label for the legend
  nbap.m$value2<-cut(nbap.m$value,breaks=c(-Inf, 0.001, 0.01, 0.05),label=c("***", "** ", "*  ")) 
  nbar.m<-cbind.data.frame(nbar.m,nbap.m$value,nbap.m$value2) # adding the p value and its cut to the first dataset of R coefficients
  names(nbar.m)[5]<-paste("valuep") # change the column names of the dataframe 
  names(nbar.m)[6]<-paste("signif.")
  nbar.m$row <- factor(nbar.m$row, levels=rev(unique(as.character(nbar.m$variable)))) # reorder the variable factor
  # Plotting the matrix correlation heatmap
  # Set options for a blank panel
  po.nopanel <-list(opts(panel.background=theme_blank(),panel.grid.minor=theme_blank(),panel.grid.major=theme_blank()))
  pa<-ggplot(nbar.m, aes(row, variable)) +
  geom_tile(aes(fill=value2),colour="white") +
  scale_fill_brewer(palette = "RdYlGn",name="Correlation")+ # RColorBrewer package
  opts(axis.text.x=theme_text(angle=-90))+
  po.nopanel
  pa # check the first plot
  # Adding the significance level stars using geom_text 
  pp<- pa +
  geom_text(aes(label=signif.),size=2,na.rm=TRUE) # you can play with the size
  # Workaround for the alpha aesthetics if it is good to represent significance level, the same workaround can be applied for size aesthetics in ggplot2 as well. Applying the alpha aesthetics to show significance is a little bit problematic, because we want the alpha to be low while the p value is high, and vice verse which can't be done without a workaround
  nbar.m$signif.<-rescale(as.numeric(nbar.m$signif.),to=c(0.1,0.9)) # I tried to use to=c(0.1,0.9) argument as you might expect, but to avoid problems with the next step of reciprocal values when dividing over one, this is needed for the alpha aesthetics as a workaround
  nbar.m$signif.<-as.factor(0.09/nbar.m$signif.) # the alpha now behaves as wanted  except for the NAs values stil show as if with three stars level, how to fix that?
  # Adding the alpha aesthetics in geom_point in a shape of squares (you can improve here)
  pp<- pa +
  geom_point(data=nbar.m,aes(alpha=signif.),shape=22,size=5,colour="darkgreen",na.rm=TRUE,legend=FALSE) # you can remove this step, the result of this step is seen in one of the layers in the above green heatmap, the shape used is 22 which is again a square but the size you can play with it accordingly
  

  我希望这可以向前迈进一步！请注意：

  • 有人建议对R ^ 2进行不同的分类或切割，好吧我们当然可以做到这一点但我们仍然希望向 Spectator 展示显着性水平，而不是用明星级别来扰乱眼睛 . 我们原则上是否可以实现这一点？
  • 有人建议用不同的方法来削减p值，好吧这可以是在没有显示3级重要性而不会使眼睛麻烦的情况下的选择 . 然后，没有水平显示显着/不显着可能更好
  • 对于alpha和size美学，你可能会更好地想出ggplot2中的上述解决方法，希望很快收到你的来信！
  • 问题尚未得到解答，等待创新的解决方案！ - 有趣的是，"corrplot"包做到了！我通过这个包得出了下面的图，PS：交叉的正方形不重要，signif = 0.05 . 但是我们怎样才能把它翻译成ggplot2，我们可以吗？！

  

  • 或者你可以做圈子并隐藏那些不重要的圈子？怎么在ggplot2中这样做？！

  

  回复于 2024-05-13T08:34:45+08:00
• 1

  为了表示估计相关系数的显着性，您可以改变着色量 - 使用 alpha 或仅填充每个图块的子集：

  # install.packages("fdrtool")
  # install.packages("data.table")
  library(ggplot2)
  library(data.table)
  
  #download dataset
  nba <- as.matrix(read.csv("http://datasets.flowingdata.com/ppg2008.csv")[-1])
  m <- ncol(nba)
  # compute corellation and p.values for all combinations of columns
  dt <- CJ(i=seq_len(m), j=seq_len(m))[i<j]
  dt[, c("p.value"):=(cor.test(nba[,i],nba[,j])$p.value), by=.(i,j)]
  dt[, c("corr"):=(cor(nba[,i],nba[,j])), by=.(i,j)]
  
  # estimate local false discovery rate
  dt[,lfdr:=fdrtool::fdrtool(p.value, statistic="pvalue")$lfdr]
  
  dt <- rbind(dt, setnames(copy(dt),c("i","j"),c("j","i")), data.table(i=seq_len(m),j=seq_len(m), corr=1, p.value=0, lfdr=0))
  
  
  #use alpha
  ggplot(dt, aes(x=i,y=j, fill=corr, alpha=1-lfdr)) + 
    geom_tile()+
    scale_fill_distiller(palette = "RdYlGn", direction=1, limits=c(-1,1),name="Correlation") +
    scale_x_continuous("variable", breaks = seq_len(m), labels = colnames(nba)) +
    scale_y_continuous("variable", breaks = seq_len(m), labels = colnames(nba), trans="reverse") +
    coord_fixed() +
    theme(axis.text.x=element_text(angle=90, vjust=0.5),
          panel.background=element_blank(),
          panel.grid.minor=element_blank(),
          panel.grid.major=element_blank(),
    )
  

  

  #use area
  ggplot(dt, aes(x=i,y=j, fill=corr,  height=sqrt(1-lfdr),  width=sqrt(1-lfdr))) + 
    geom_tile()+
    scale_fill_distiller(palette = "RdYlGn", direction=1, limits=c(-1,1),name="Correlation") +
    scale_color_distiller(palette = "RdYlGn", direction=1, limits=c(-1,1),name="Correlation") +
    scale_x_continuous("variable", breaks = seq_len(m), labels = colnames(nba)) +
    scale_y_continuous("variable", breaks = seq_len(m), labels = colnames(nba), trans="reverse") +
    coord_fixed() +
    theme(axis.text.x=element_text(angle=90, vjust=0.5),
          panel.background=element_blank(),
          panel.grid.minor=element_blank(),
          panel.grid.major=element_blank(),
    )
  

  

  这里的关键是p.values的缩放：为了获得仅在相关区域中显示大变化的易于理解的值，我使用由 fdrtools 提供的本地错误发现（lfdr）的上限估计 . 即，图块的alpha值可能小于或等于该相关的概率与0不同 .

  回复于 2024-05-13T08:34:45+08:00