从包含许多二元分类列的 data.table 中观察的流行率估计

Question

我正在通过蛮力从我的原始 data.table 进行流行率估计，我需要提高效率。你能帮忙吗？

我的 data.table 每行包含一个加权观察值。有许多列充当二进制虚拟变量，指示特定观察是否属于许多可能分类中的一个或多个。 （例如，一个故事可以是“惊人的”、“无聊的”或“迷人的”，或三者的任意组合。）

必须有一个data.table 的方式来替换我的forloop。我还怀疑我可能不需要生成queries 集。我很感激以全新的眼光看待这个问题。

library(data.table)

set.seed(42)
# I have many weighted observations that can be labeled as belonging to one of many categories
# in this example, I simulate 10 observations and only 3 categories
dt = data.table(
        weight = runif( n = 10 , min = 0, max = 1 ),
        a = sample( x = c(0,1) , size = 10 , replace = TRUE ),
        b = sample( x = c(0,1) , size = 10 , replace = TRUE ),
        c = sample( x = c(0,1) , size = 10 , replace = TRUE )
)

# Generate all combinations of categories
queries = as.data.table( expand.grid( rep( list(0:1) , length(names(dt))-1 ) ) )
names(queries) = names(dt)[ 2:length(names(dt)) ] # rename Var1, Var2, Var3 to a, b, c

# Brute force through each possible combination to calculate prevalence
prevalence = rep( NA, nrow(queries) )
for( q in 1:nrow(queries) ){
    prevalence[q] = dt[ a == queries[q, a] & b == queries[q, b] & c == queries[q, c] , sum(weight) ] / dt[ , sum(weight) ]
}

results = copy(queries)
results$prevalence = prevalence

results

输出是：

#   a b c prevalence
#1: 0 0 0 0.09771385
#2: 1 0 0 0.10105192
#3: 0 1 0 0.36229784
#4: 1 1 0 0.00000000
#5: 0 0 1 0.00000000
#6: 1 0 1 0.05993197
#7: 0 1 1 0.00000000
#8: 1 1 1 0.37900443

更新：原始问题有 42 个模拟观察结果，数据涵盖了每个可能的类别组合（a、b、c）。该问题已修改为仅包含 10 个模拟观察结果，因此会有没有观察结果的组合（且流行率为零）。

Answer 1

更新答案

方法一：

1. 使用CJ创建a,b,c的完整组合，然后加入dt（如@TMo's answer）
2. 将每组的weight相加，然后除以totoal_weight
3. NA 的出现是理性的。如果需要，您也可以使用nafill 函数填充0。

total_weight = sum(dt$weight)
dt[CJ(a, b, c, unique = TRUE),  
   on = .(a, b, c)][, 
                    .( prevalence = sum(weight)/total_weight), 
                    by = .(a,b,c)]

#      a     b     c prevalence
#   <num> <num> <num>      <num>
#1:     0     0     0 0.09771385
#2:     0     0     1         NA
#3:     0     1     0 0.36229784
#4:     0     1     1         NA
#5:     1     0     0 0.10105192
#6:     1     0     1 0.05993197
#7:     1     1     0         NA
#8:     1     1     1 0.37900443

方法二：

dt2 = dt[,.( prevalence = sum(weight) / total_weight ), by = .(a,b,c)]
dt2[queries, on = .(a,b,c)]
# or `queries[, prevalence := fcoalesce(dt2[queries, prevalence])]`
#       a     b     c prevalence
#   <int> <int> <int>      <num>
#1:     0     0     0 0.09771385
#2:     1     0     0 0.10105192
#3:     0     1     0 0.36229784
#4:     1     1     0         NA
#5:     0     0     1         NA
#6:     1     0     1 0.05993197
#7:     0     1     1         NA
#8:     1     1     1 0.37900443

原答案

可以按组计算

dt[,.( prevalence = sum(weight) / dt[,sum(weight)] ), by = .(a,b,c)]

• 每个组对应您的类别
• 将每组的weight相加，然后除以总权重

Answer 2

这里有一些解决方案（在这两种情况下，您都可以将keyby 参数替换为by）

如果您的数据集 (dt) 已经包含不同类别的所有可能组合，那么您可以这样做（如 @Peace Wang solution）

dt[, .(prevalence = sum(weight)/sum(dt$weight)), keyby=.(a, b, c)]

#        a     b     c prevalence
# 1:     0     0     0 0.10876301
# 2:     0     0     1 0.02135357
# 3:     0     1     0 0.03775363
# 4:     0     1     1 0.12806864
# 5:     1     0     0 0.18204696
# 6:     1     0     1 0.15197811
# 7:     1     1     0 0.25629705
# 8:     1     1     1 0.11373903

相反，如果数据集不包含不同类别的所有可能组合，那么您可以如下解决它（CJ(a, b, c, unique=TRUE) 计算所有组合并删除重复项）

dt[CJ(a, b, c, unique=TRUE), .(prevalence = sum(weight)/sum(dt$weight)), keyby=.(a, b, c), on=.(a, b, c)]

#        a     b     c prevalence
# 1:     0     0     0 0.10876301
# 2:     0     0     1 0.02135357
# 3:     0     1     0 0.03775363
# 4:     0     1     1 0.12806864
# 5:     1     0     0 0.18204696
# 6:     1     0     1 0.15197811
# 7:     1     1     0 0.25629705
# 8:     1     1     1 0.11373903