R语言绘制一次和二次相关性热图

  在数据探索的过程中，我们往往会对数据与数据的相关性进行分析，例如我们常用的corrplot包，或者psych包中的corr.test函数，对两两变量间的相关性进行分析。我们常常会看到这样的相关性热图：

  但有时变量间的关系并非线性，例如：我们熟知的微生物物种多样性与土壤pH之间的关系，因此在数据探索过程中，仅使用线性相关会导致我们遗漏一些重要信息。二次关系也需要纳入我们考虑的范围，因此本篇博客将对两两变量的二次关系的R代码及可视化进行介绍，助力数据探索。

  我们尝试在对数据进行线性相关分析的同时，加入二次关系的分析，数据的可视化显示的结果大致如下：其中菱形表示线性相关显著，圆圈表示线性相关不显著，而二次关系显著。

可视化代码遵循如下原则：

• 若线性相关性显著，无论二次关系是否显著，均使用线性相关的R值和P值；（使用最简单的模型即可）
• 若线性相关性不显著，则查看二次关系是否显著，若显著则显示二次关系的R值和P值；若不显著则均不显示。

1.计算变量的两两二次关系函数

qcorr <- function(x){
  pm <- matrix(rep(0,ncol(x)^2),nrow = ncol(x))
  rm <- matrix(rep(0,ncol(x)^2),nrow = ncol(x))
  rownames(pm) <- colnames(x)
  colnames(pm) <- colnames(x)
  rownames(rm) <- colnames(x)
  colnames(rm) <- colnames(x)
  for(i in 1:length(x)){
    for(j in 1:length(x)){
      if (i != j){
        x1 <- x[,i]
        x2 <- x[,j]
        s <- summary(lm(x1~x2+I(x2^2)))
        pm[i,j] <- pf(s$fstatistic[1],s$fstatistic[2],s$fstatistic[3],lower.tail=F)
        rm[i,j] <- sqrt(s$r.squared)
      }
    }
  }
  return(list(p=pm,r=rm))
}

2.与线性相关性进行整合

env_corr2 <- qcorr(env_div)$r
env_corp2 <- qcorr(env_div)$p
env_corr1 <- corr.test(env_div,adjust = "none")$r
env_corp1 <- corr.test(env_div,adjust = "none")$p

corr2 <- melt(env_corr2, na.rm = TRUE)
corp2 <- melt(env_corp2, na.rm = TRUE)
corr1 <- melt(env_corr1, na.rm = TRUE)
corp1 <- melt(env_corp1, na.rm = TRUE)
type <- ifelse(corp1$value>0.05 & corp2$value < 0.05,"Q","L")
corp1[corp1$value>0.05 & corp2$value < 0.05,] <- corp2[corp1$value >0.05 & corp2$value < 0.05,]
corr1[corp1$value>0.05 & corp2$value < 0.05,] <- corr2[corp1$value >0.05 & corp2$value < 0.05,]
corr <- corr1
corr$p_value <- corp1$value
corr$type <- type
corr$abs_corr <- abs(corr$value)*10
corr <- corr[corr$p_value<0.05 & corr$value < 1,]

3.对整合的结果进行可视化

outline.color = "gray"
colors = c("blue", "white", "red")
shape <- cols <- c("Q" = 21, "L" = 23)
library(ggplot2)
p <- ggplot(data = corr, mapping = aes_string(x = "Var1", y = "Var2"))
p <- p + geom_point(color = outline.color,aes_string(size = "abs_corr", fill = "value",shape = "type")) +
         scale_size(range = c(4,10)) +
         scale_shape_manual(values = shape) +
         guides(size = "none")
p <- p + scale_fill_gradient2(low = colors[1], high = colors[3], 
                              mid = colors[2], midpoint = 0, limit = c(-1, 1), space = "Lab")
p <- p + ggplot2::theme_minimal()
p

这里使用的是ggplot2的包，方便对可视化的结果进一步自定义。
输出结果如下：

测试数据及代码可参考如下连接：https://mbd.pub/o/bread/ZZ6WmZhr