从方阵用 R 创建热图

Creating heat map with R from a square matrix

我有一个 gzip 压缩文件 file.gz,其中包含 4,726,276 行,其中第一行和最后五行如下所示:

 FID1 IID1 FID2 IID2 RT    EZ      Z0      Z1      Z2  PI_HAT PHE       DST     PPC   RATIO
 CAN  -1 CAN   1 OT     0  1.0000  0.0000  0.0000  0.0000  -1  0.745118  0.1111  1.5526
 CAN  -1 CAN   2 OT     0  0.8761  0.1239  0.0000  0.0619  -1  0.752607  0.0648  1.4615
 CAN  -1 CAN   3 OT     0  0.8810  0.1190  0.0000  0.0595  -1  0.753934  0.3058  1.7941
 CAN  -1 CAN   4 OT     0  0.8911  0.1089  0.0000  0.0545  -1  0.751706  0.8031  2.4138

 WAN   2 WAN   4 OT     0  0.8410  0.0000  0.1590  0.1590  -1  0.787251  0.0840  1.5000
 WAN   2 WAN   5 OT     0  0.8606  0.0000  0.1394  0.1394  -1  0.784882  0.7671  2.3571
 WAN   3 WAN   4 OT     0  0.8306  0.0000  0.1694  0.1694  -1  0.790142  0.0392  1.3846
 WAN   3 WAN   5 OT     0  0.7960  0.0364  0.1676  0.1858  -1  0.795924  0.5000  2.0000
 WAN   4 WAN   5 OT     0  0.8227  0.0090  0.1683  0.1728  -1  0.793460  0.5577  2.0645

x 值是第 1 + 2 列。y 值是第 3 + 4 列。z 值是第 10 列。输入文件中不存在沿对角线的值。它们最好是 1,但 0 也可以。

如何根据此类数据创建热图?

这是一个 3x3 矩阵的简单示例:

FID1 IID1 FID2 IID2 PI_HAT
A 1 B 1 0.1
A 1 B 2 0.2
B 1 B 2 0.3

您的问题似乎分为两部分:

  1. 如何处理 R 中的数据(在本例中来自 gzip 压缩存档)
  2. 如何制作热图

乍一看,您似乎在暗示数据量很大——有很多关于如何在 R 中使用大数据的资源 (here's one)——但是基于评论 我认为数据大小实际上不是问题。如果是,那么您的选择将部分取决于您的硬件资源以及您采样数据的意愿(我强烈推荐),而不是使用 500 万行中的每一行。 Central Limit Theorem 是你的朋友。

您可以像这样读入 gzip 数据:

data <- read.table(gzfile("file.gz"),header=T, sep="\t", stringsAsFactors=F)

由于您没有提供压缩档案,我已经复制了您的示例数据并在下面的代码中从我的剪贴板中读取了它。我将向您展示如何根据这些数据构建热图;要从 gzip 导入和处理大数据,请查看上面提供的 link。

require(stats)
require(fields)
require(akima)

a <- read.table(con <- file("clipboard"), header = T)

a$x1 <- as.numeric(a[,1])
a$x2 <- as.numeric(a[,2])
a$y1 <- as.numeric(a[,3])
a$y2 <- as.numeric(a[,4])
x <- as.matrix(cbind(a$x1, a$x2))
y <- as.matrix(cbind(a$y1, a$y2))
z <- as.matrix(a[, 10])

s = smooth.2d(z, x=cbind(x,y), theta=0.5)
image.plot(s)

这是一种ggplot2方法。 450 万行在 R 中应该不是问题。

df <- read.table(text='FID1 IID1 FID2 IID2 RT    EZ      Z0      Z1      Z2  PI_HAT PHE       DST     PPC   RATIO
 CAN  -1 CAN   1 OT     0  1.0000  0.0000  0.0000  0.0000  -1  0.745118  0.1111  1.5526
 CAN  -1 CAN   2 OT     0  0.8761  0.1239  0.0000  0.0619  -1  0.752607  0.0648  1.4615
 CAN  -1 CAN   3 OT     0  0.8810  0.1190  0.0000  0.0595  -1  0.753934  0.3058  1.7941
 CAN  -1 CAN   4 OT     0  0.8911  0.1089  0.0000  0.0545  -1  0.751706  0.8031  2.4138
 CAN  -1 CAN   4 OT     0  0.8911  0.1089  0.0000  0.0545  -1  0.751706  0.8031  2.4138
 CAN  -1 CAN   4 OT     0  0.8911  0.1089  0.0000  0.0545  -1  0.751706  0.8031  2.4138
 WAN   3 WAN   4 OT     0  0.8306  0.0000  0.1694  0.1694  -1  0.790142  0.0392  1.3846
 WAN   3 WAN   4 OT     0  0.8306  0.0000  0.1694  0.1694  -1  0.790142  0.0392  1.3846
 WAN   2 WAN   4 OT     0  0.8410  0.0000  0.1590  0.1590  -1  0.787251  0.0840  1.5000
 WAN   2 WAN   5 OT     0  0.8606  0.0000  0.1394  0.1394  -1  0.784882  0.7671  2.3571
 WAN   3 WAN   4 OT     0  0.8306  0.0000  0.1694  0.1694  -1  0.790142  0.0392  1.3846
 WAN   3 WAN   5 OT     0  0.7960  0.0364  0.1676  0.1858  -1  0.795924  0.5000  2.0000
 WAN   4 WAN   5 OT     0  0.8227  0.0090  0.1683  0.1728  -1  0.793460  0.5577  2.0645', header=T)

我在您的输出中添加了几行,以使其在热图中更合理。之前没有重叠:

#create your variables by merging columns 1+2 and 3+4
a <- mapply(paste,df[[1]], df[[2]])
b <- mapply(paste,df[[3]], df[[4]])
#combine in a data.frame
df2 <- data.frame(a,b)


library(dplyr)
#aggregate because you will need aggregated rows for this to work
#this should only take a few seconds for 4.5m rows
df3 <-
df2 %>%
  group_by(a,b) %>%
  summarize(total=n())

#plot with ggplot2
library(ggplot2)
ggplot(df3, aes(x=a,y=b,fill=total)) + geom_tile()

输出: