考察一個群體的多個表型或者一個表型的多個重復形娇,我們想展示其分布和他們之間的相關關系可以使用柱狀圖和散點圖(如下圖所示)抠忘。
這幅圖主要有兩部分組成濒生,一個是對角線上的柱狀圖,使用柱狀圖展示了每一個表型重復的分布砂竖;另一個就是對角線下面的散點圖,用散點圖展示兩兩之間的相關關系鹃答,并且用不同顏色表示點的密度乎澄,在上面標注其相關性。下面我們將使用R語言完成這幅圖测摔。
對于這幅圖我們可以先分別繪制其中每一個部分置济,然后使用圖片組合、拼接函數(shù)進行整合:
分圖繪制
首先導入數(shù)據(jù)锋八,數(shù)據(jù)格式如下浙于,每一行代表一個樣本,每一列代表一個重復:
>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()
函數(shù)繪制,使用cor()
函數(shù)計算兩個重復之間的相關系數(shù)樊销,并將其放在圖片標題位置整慎,并使用ggtext
包的element_markdown()
函數(shù)設置標題的主題,同時使用cowplot
包的theme_half_open()
函數(shù)設置整體主題围苫。
>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"))
組合圖片
使用customLayout包進行圖片組合裤园,這個包可以對base繪圖和ggplot2繪圖進行整合,而且比較靈活剂府。首先需要lay_new()函數(shù)創(chuàng)建一個拼接畫布拧揽,然后使用lay_grid()函數(shù)組合各個圖片。因為總共有5個重復腺占,因此需要一個5×5的畫圖淤袜,如下圖所示,各個分圖從左上角開始往下排列走”之“字形排列衰伯。
>lay <- lay_new(mat = matrix(1:25, nrow = 5), widths = rep(1, 1), heights = rep(1, 1))
>lay_show(lay)
現(xiàn)在出現(xiàn)了一個問題铡羡,我們并沒有在對角線上方安排圖片,而lay_new()產(chǎn)生的是一個矩形排列畫布意鲸,因此我們需要在右上角填充空白圖片烦周,并將空白圖和柱狀圖、散點密度圖整合怎顾。
>p <- ggplot() + theme_nothing()
>lay_grid(list(p1, p2, p3, ...), lay)
整理以上過程
在一個5×5的組合中我們總共需要繪制25個分圖读慎,其中有多次重復的過程,并且最終圖片是矩形有規(guī)律分布槐雾,因此為了減少代碼長度我們可以使用循環(huán)來處理每個分圖夭委。根據(jù)lay_new()的組合形式可以設置兩層循環(huán)分別處理行和列,并且因為組合圖是從左上角開始向下排布募强,因此外層循環(huán)用來處理行株灸,內層分布處理列崇摄。最后一點就是可以把這一系列代碼寫成一個function,方便以后使用蚂且。
最終代碼如下所示:
library(tidyverse)
library(ggpointdensity)
library(cowplot)
library(ggtext)
library(customLayout)
## 定義pheCorDist函數(shù)
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)
}
# 導入數(shù)據(jù)并繪圖
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()