R语言R学习与可视化走进转录组

R绘图:群体表型分布和相关关系

2022-06-10  本文已影响0人  菜刀_ebc3

考察一个群体的多个表型或者一个表型的多个重复,我们想展示其分布和他们之间的相关关系可以使用柱状图和散点图(如下图所示)。


test.png

这幅图主要有两部分组成,一个是对角线上的柱状图,使用柱状图展示了每一个表型重复的分布;另一个就是对角线下面的散点图,用散点图展示两两之间的相关关系,并且用不同颜色表示点的密度,在上面标注其相关性。下面我们将使用R语言完成这幅图。

对于这幅图我们可以先分别绘制其中每一个部分,然后使用图片组合、拼接函数进行整合:

分图绘制

首先导入数据,数据格式如下,每一行代表一个样本,每一列代表一个重复:

>data <- read.table("./data.txt", header = T, row.names = 1, sep = "\t")
>data
     pheno16rep1 pheno16rep2 pheno17 pheno18rep1 pheno18rep2
s1          28.4        24.9  27.740    27.72500    29.30000
s2          26.6        25.3      NA          NA          NA
s3          27.8        27.0  24.660          NA    27.97500
s4          25.5        26.9  22.680    29.27500    25.95000
s5          26.5        28.7  24.760    31.97500    27.52000
s6          ....        ....  .......   ........    ........

使用ggplot2扩展包绘制每一个分图。柱状图使用geom_histogram()绘制,散点图使用ggpointdensity包的geom_pointdensity()函数绘制,使用cor()函数计算两个重复之间的相关系数,并将其放在图片标题位置,并使用ggtext包的element_markdown()函数设置标题的主题,同时使用cowplot包的theme_half_open()函数设置整体主题。

>library(tidyverse)
>library(ggpointdensity)
>library(cowplot)
>library(ggtext)
>
># 柱状图
>p1 <- ggplot(data, aes(x = pheno16rep1)) + 
>        geom_histogram(binwidth = 1) +
>        labs(x = NULL, y = NULL, title = cn[i]) +
>        theme_half_open() +
>        theme(plot.title = element_text(hjust = 0.5))
># 散点密度图
>p2 <- ggplot(data, aes(x = pheno16rep2, y = pheno16rep1)) +
>        geom_pointdensity() +
>        scale_color_continuous(type = "viridis") +                # 设置点密度颜色梯度
>        labs(x = NULL, y = NULL, title = paste("*R*: ", round(cor(df$col1, df$col2, use = "na.or.complete"), 2), sep = "")) +
>        theme_half_open() +
>        theme(legend.position = "NA", 
>           plot.title = element_markdown(hjust = 0.5,
>                           face = "plain"))
Rplot03.png

组合图片

使用customLayout包进行图片组合,这个包可以对base绘图和ggplot2绘图进行整合,而且比较灵活。首先需要lay_new()函数创建一个拼接画布,然后使用lay_grid()函数组合各个图片。因为总共有5个重复,因此需要一个5×5的画图,如下图所示,各个分图从左上角开始往下排列走”之“字形排列。

>lay <- lay_new(mat = matrix(1:25, nrow = 5), widths = rep(1, 1), heights = rep(1, 1))
>lay_show(lay)
Rplot02.png

现在出现了一个问题,我们并没有在对角线上方安排图片,而lay_new()产生的是一个矩形排列画布,因此我们需要在右上角填充空白图片,并将空白图和柱状图、散点密度图整合。

>p <- ggplot() + theme_nothing()
>lay_grid(list(p1, p2, p3, ...), lay)

整理以上过程

在一个5×5的组合中我们总共需要绘制25个分图,其中有多次重复的过程,并且最终图片是矩形有规律分布,因此为了减少代码长度我们可以使用循环来处理每个分图。根据lay_new()的组合形式可以设置两层循环分别处理行和列,并且因为组合图是从左上角开始向下排布,因此外层循环用来处理行,内层分布处理列。最后一点就是可以把这一系列代码写成一个function,方便以后使用。

最终代码如下所示:

library(tidyverse)
library(ggpointdensity)
library(cowplot)
library(ggtext)
library(customLayout)

## 定义pheCorDist函数
pheCorDist <- function(data) {
  #
  n <- ncol(data)
  cn <- colnames(data)
  Pall <- list()
  index <- 1
  #
  for (j in 1:n) {
    for (i in 1:n) {
      df <- data[, c(i,j)]
      colnames(df) <- c("col1", "col2")
      if (i == j) {
        p <- ggplot(df, aes(x = col1)) + 
          geom_histogram(binwidth = 1) +
          labs(x = NULL, y = NULL, title = cn[i]) +
          theme_half_open() +
          theme(plot.title = element_text(hjust = 0.5))
      } else if (i > j) {
        p <- ggplot(df, aes(x = col2, y = col1)) +
          geom_pointdensity() +
          scale_color_continuous(type = "viridis") +
          labs(x = NULL, y = NULL, title = paste("*R*: ", round(cor(df$col1, df$col2, use = "na.or.complete"), 2), sep = "")) +
          theme_half_open() +
          theme(legend.position = "NA", 
                plot.title = element_markdown(hjust = 0.5,
                                              face = "plain"))
      } else if(i < j) {
        p <- ggplot() + theme_nothing()
      }
      Pall[index][[1]] <- p
      index = index + 1
    }
  }
  lay <- lay_new(mat = matrix(1:n^2, nrow = n), widths = rep(1, n), heights = rep(1, n))
  lay_grid(Pall, lay)
}

# 导入数据并绘图
data <- read.table("./data.txt", header = T, row.names = 1, sep = "\t")
png(filename = "test.png", width = 10, height = 8, units = "in", res = 500)
pheCorDist(data)
dev.off()
上一篇下一篇

猜你喜欢

热点阅读