logo标签列表

按类别计算列出现次数

rcountdata.table
3
我正在尝试计算每个SNP名称中iets列中“Opp”出现的次数（最终我想将“Opp”的出现次数除以df$MM）。
library(data.table)
df <- structure(list(SNP = structure(c(1L, 1L, 1L, 2L, 1L), .Label = c("rs80932150", "rs000001"), class = "factor"), FID = c(116601888L, 116621563L, 117253533L, 118635095L, 118943247L), IID = c(116601888L, 116621563L, 117253533L, 118635095L, 118943247L), NEW = structure(c(16L, 14L, 16L, 14L, 14L), .Label = c("A/A", "A/C", "A/G", "A/T", "C/A", "C/C", "C/G", "C/T", "G/A", "G/C", "G/G", "G/T", "T/A", "T/C", "T/G", "T/T"), class = "factor"), OLD = structure(c(6L, 6L, 6L, 6L, 6L), .Label = c("A/A", "A/C", "A/G", "A/T", "C/A", "C/C", "C/G", "C/T", "G/A", "G/C", "G/G", "G/T", "T/A", "T/C", "T/G", "T/T"), class = "factor"), count = c(1L, 1L, 1L, 1L, 1L), MM = c(4L, 4L, 4L, 1L, 4L), iets = c("Opp", "Het", "Opp", "Het", "Het")), .Names = c("SNP", "FID", "IID", "NEW", "OLD", "count", "MM", "iets"), class = "data.frame", row.names = c(NA, -5L))
setDT(df)

#         SNP       FID       IID NEW OLD count MM iets
#1 rs80932150 116601888 116601888 T/T C/C     1  4  Opp
#2 rs80932150 116621563 116621563 T/C C/C     1  4  Het
#3 rs80932150 117253533 117253533 T/T C/C     1  4  Opp
#4   rs000001 118635095 118635095 T/C C/C     1  1  Het
#5 rs80932150 118943247 118943247 T/C C/C     1  4  Het

我期望的结果如下：

df
#          SNP       FID       IID NEW OLD count MM iets oppcount percentage
#1: rs80932150 116601888 116601888 T/T C/C     1  4  Opp      2        0.5
#2: rs80932150 116621563 116621563 T/C C/C     1  4  Het      2        0.5
#3: rs80932150 117253533 117253533 T/T C/C     1  4  Opp      2        0.5
#4:   rs000001 118635095 118635095 T/C C/C     1  1  Het      0        0.0
#5: rs80932150 118943247 118943247 T/C C/C     1  4  Het      2        0.5

我一直在尝试类似的事情，但是我似乎无法弄清楚如何将出现值分配给我的oppcount/percentage列。
首先，我需要计算每个SNP中"Opp"的数量，然后将其除以MM。

as.character((sum(df$iets == "Opp")/(df[,.N, by = df$SNP][[2]])))
#[1] "0.5" "2"

我该如何计算每个SNP（类别）中“Opp”出现的次数？
3个回答
4
您可以使用:=操作符通过引用更新您的data.table。 例如：
df[, `:=` (oppcount = sum(iets=='Opp'), percentage = sum(iets=='Opp')/.N), by = SNP]

你需要知道的是:

> df
          SNP       FID       IID NEW OLD count MM iets oppcount percentage
1: rs80932150 116601888 116601888 T/T C/C     1  4  Opp        2        0.5
2: rs80932150 116621563 116621563 T/C C/C     1  4  Het        2        0.5
3: rs80932150 117253533 117253533 T/T C/C     1  4  Opp        2        0.5
4:   rs000001 118635095 118635095 T/C C/C     1  1  Het        0        0.0
5: rs80932150 118943247 118943247 T/C C/C     1  4  Het        2        0.5

或者，根据评论中@Frank的建议，您还可以使用以下两个选项之一：

# method 1
df[, c('oppcount', 'percentage') := {s = sum(iets=='Opp'); .(s, s/.N)}, by = SNP]
# method 2
df[df[, {s = sum(iets=='Opp'); .(oppcount = s, percentage = s/.N)}, by = SNP], on = 'SNP']

一个基础的R替代方案：
transform(df,
          oppcount = ave(iets, SNP, FUN = function(x) sum(x=='Opp')),
          percentage = ave(iets, SNP, FUN = function(x) sum(x=='Opp')/length(x)))

一个正确的dplyr替代方案为:

library(dplyr)
df %>% 
  group_by(SNP) %>% 
  mutate(oppcount = sum(iets=='Opp'),
         percentage = oppcount/n())
0
rs8.oppcount<-length(iets[iets=='Opp' & SNP=='rs80932150'])
rs0.oppcount<-length(iets[iets=='Opp' & SNP=='rs000001'])

这将保存 Opp 在 snp 类别中的出现次数！

编辑：

df1<-group_by(df, df$SNP)
df2<-summarise(df1, oppcount = length(iets[iets=='Opp']))
df1<-merge(df1, df2, by = 'SNP')

这个可以吗？

这是第一步，但我的实际数据集包含大约10万行，总共有约20K个不同的SNP。我正在寻找一种无论SNP名称如何都能工作的解决方案。 - Bas
我希望按照SNP组重新排序不会有影响。 - Daniel Winkler
很遗憾，代码出了问题。> df2<-summarise(df1, oppcount = length(iets[iets=='Opp'])) > df1<-merge(df1, df2, by = 'SNP') Error in merge.data.table(df1, df2, by = "SNP") : Elements listed in by must be valid column names in x and y. 我无法让它正常工作，因为df2有一个名为df2$df$SNP的列。 - Bas
重新排序SNP组并不重要，我可以简单地使用match()将它们放回原位！ - Bas
顺便提一下，计算满足条件的实例数量的标准方法是 sum(cond(x)) 而不是 length(x[ cond(x) ])，这在被接受的答案中可以看到。我认为 Dplyr 有额外的方便函数来处理这个问题。 - Frank
0

使用dplyr怎么样？

library('dplyr')
df %>% group_by(iets, SNP) %>% summarize(count=sum(count)) %>% filter(iets=='Opp')
1
这只返回一个包含示例的行，不包括仅为“Het”的rs id。 - Bas
网页内容由stack overflow 提供, 点击上面的
可以查看英文原文，
原文链接