R中基因型的设计矩阵
Design matrix for genotypes in R
我正在寻找一种有效的方法来为 R 中的基因型创建“参数化”设计矩阵。
我有一个大文件(大约 3 GB),里面有动物和它们的基因型。示例数据如下所示:
snp id a1 a2 code
snp1 an1 A A 0
snp1 an2 A B 1
snp1 an3 B B -1
snp2 an1 A B 1
snp2 an2 A A 0
snp2 an3 B B -1
snp是snp的名称(每只动物都有各自的snp),id是动物的id(每只动物都有唯一的id),a1是等位基因1,a2是等位基因2,代码表示基于等位基因的基因型,如果动物有两个A它的代码是0,如果动物有AB,它的代码是-1,如果它是BB的代码是1.
现在我需要基于该设计矩阵创建,它在行中有动物(数据中的 id 列),在 SNP 列中(数据中的 snp 列)和“单元格”(在交集处)行和列)我需要代码列中的值。所以最后,它应该是这样的:
an1 0 1
an2 1 0
an3 -1 -1
我知道在效率和速度方面 R 有限制,但我仍然需要 R 中可以获得的最快解决方案。
通常 data.table 包在这些类型的情况下性能非常好。示例如下:
library(data.table)
#> Warning: package 'data.table' was built under R version 4.1.1
df <- fread(text = "snp id a1 a2 code
snp1 an1 A A 0
snp1 an2 A B 1
snp1 an3 B B -1
snp2 an1 A B 1
snp2 an2 A A 0
snp2 an3 B B -1")
dcast(df, id ~ snp, value.var = "code")
#> id snp1 snp2
#> 1: an1 0 1
#> 2: an2 1 0
#> 3: an3 -1 -1
由 reprex package (v2.0.1)
于 2021-10-13 创建
如果您需要将输出作为矩阵,您可以使用:
cast <- dcast(df, id ~ snp, value.var = "code")
mat <- as.matrix(cast[, -"id"])
rownames(mat) <- cast$id
mat
#> snp1 snp2
#> an1 0 1
#> an2 1 0
#> an3 -1 -1
对于大约 3Gb 的文件,您可能希望 运行 持续大约 10 秒:
library(data.table)
#> Warning: package 'data.table' was built under R version 4.1.1
# Setting up larger data
df <- expand.grid(
snp = paste0("snp", 1:10000),
id = paste0("an", 1:10000)
)
df$a1 <- sample(c("A", "B"), nrow(df), replace = TRUE)
df$a2 <- sample(c("A", "B"), nrow(df), replace = TRUE)
df$code <- with(df, dplyr::case_when(
a1 == "A" & a2 == "A" ~ 0,
a1 == "B" & a2 == "B" ~ -1,
TRUE ~ 1
))
setDT(df)
# How big is this data?
format(object.size(df), "Gb")
#> [1] "3 Gb"
# How fast does the function run?
bench::mark(
dcast(df, id ~ snp, value.var = "code")
)
#> Warning: Some expressions had a GC in every iteration; so filtering is disabled.
#> # A tibble: 1 x 6
#> expression min median `itr/sec` mem_alloc
#> <bch:expr> <bch:tm> <bch:tm> <dbl> <bch:byt>
#> 1 dcast(df, id ~ snp, value.var = "code") 9.32s 9.32s 0.107 6.71GB
#> # ... with 1 more variable: gc/sec <dbl>
由 reprex package (v2.0.1)
于 2021-10-13 创建
我正在寻找一种有效的方法来为 R 中的基因型创建“参数化”设计矩阵。 我有一个大文件(大约 3 GB),里面有动物和它们的基因型。示例数据如下所示:
snp id a1 a2 code
snp1 an1 A A 0
snp1 an2 A B 1
snp1 an3 B B -1
snp2 an1 A B 1
snp2 an2 A A 0
snp2 an3 B B -1
snp是snp的名称(每只动物都有各自的snp),id是动物的id(每只动物都有唯一的id),a1是等位基因1,a2是等位基因2,代码表示基于等位基因的基因型,如果动物有两个A它的代码是0,如果动物有AB,它的代码是-1,如果它是BB的代码是1.
现在我需要基于该设计矩阵创建,它在行中有动物(数据中的 id 列),在 SNP 列中(数据中的 snp 列)和“单元格”(在交集处)行和列)我需要代码列中的值。所以最后,它应该是这样的:
an1 0 1
an2 1 0
an3 -1 -1
我知道在效率和速度方面 R 有限制,但我仍然需要 R 中可以获得的最快解决方案。
通常 data.table 包在这些类型的情况下性能非常好。示例如下:
library(data.table)
#> Warning: package 'data.table' was built under R version 4.1.1
df <- fread(text = "snp id a1 a2 code
snp1 an1 A A 0
snp1 an2 A B 1
snp1 an3 B B -1
snp2 an1 A B 1
snp2 an2 A A 0
snp2 an3 B B -1")
dcast(df, id ~ snp, value.var = "code")
#> id snp1 snp2
#> 1: an1 0 1
#> 2: an2 1 0
#> 3: an3 -1 -1
由 reprex package (v2.0.1)
于 2021-10-13 创建如果您需要将输出作为矩阵,您可以使用:
cast <- dcast(df, id ~ snp, value.var = "code")
mat <- as.matrix(cast[, -"id"])
rownames(mat) <- cast$id
mat
#> snp1 snp2
#> an1 0 1
#> an2 1 0
#> an3 -1 -1
对于大约 3Gb 的文件,您可能希望 运行 持续大约 10 秒:
library(data.table)
#> Warning: package 'data.table' was built under R version 4.1.1
# Setting up larger data
df <- expand.grid(
snp = paste0("snp", 1:10000),
id = paste0("an", 1:10000)
)
df$a1 <- sample(c("A", "B"), nrow(df), replace = TRUE)
df$a2 <- sample(c("A", "B"), nrow(df), replace = TRUE)
df$code <- with(df, dplyr::case_when(
a1 == "A" & a2 == "A" ~ 0,
a1 == "B" & a2 == "B" ~ -1,
TRUE ~ 1
))
setDT(df)
# How big is this data?
format(object.size(df), "Gb")
#> [1] "3 Gb"
# How fast does the function run?
bench::mark(
dcast(df, id ~ snp, value.var = "code")
)
#> Warning: Some expressions had a GC in every iteration; so filtering is disabled.
#> # A tibble: 1 x 6
#> expression min median `itr/sec` mem_alloc
#> <bch:expr> <bch:tm> <bch:tm> <dbl> <bch:byt>
#> 1 dcast(df, id ~ snp, value.var = "code") 9.32s 9.32s 0.107 6.71GB
#> # ... with 1 more variable: gc/sec <dbl>
由 reprex package (v2.0.1)
于 2021-10-13 创建