代码之家 › 专栏 › 技术社区 › Mauricio

在R的大数据集中按行生成多项式随机变量需要加快速度

multinomial purrr dplyr random r

Mauricio · 技术社区 · 7 年前

我有如下结构的数据,我正在尝试实现一个函数,该函数使用随机多项式生成器为每一行赋值 <Educ_W1:Educ_W5> 作为概率向量(它们逐行相加)。因此,对于每一行,新变量的值将介于1和5之间。我自己已经设法实现了它,但我试图找到一种更快的方法来实现它,因为在当前版本中它需要太长的时间(几天。数据包含100多万个观察结果)。

| IDhh|Year |Educ_W |Educ_H | Educ_W1| Educ_W2| Educ_W3| Educ_W4| Educ_W5| |----:|:----|:------|:------|---------:|---------:|---------:|---------:|---------:| | 1|1975 |2 |2 | 0.1645188| 0.3362659| 0.3940354| 0.0831637| 0.0220162| | 2|1975 |2 |2 | 0.1645188| 0.3362659| 0.3940354| 0.0831637| 0.0220162| | 5|1975 |2 |1 | 0.5103815| 0.2092249| 0.2285570| 0.0392398| 0.0125968| | 6|1975 |3 |3 | 0.0811203| 0.1535407| 0.5528233| 0.1486548| 0.0638609| | 8|1975 |1 |1 | 0.5103815| 0.2092249| 0.2285570| 0.0392398| 0.0125968| | 10|1975 |3 |2 | 0.1645188| 0.3362659| 0.3940354| 0.0831637| 0.0220162|

目前,我正在以以下方式实现该函数,但我需要非常长的时间。这里,变量“IDhh”唯一标识每一行。功能 rMultinom 来自R包 Hmisc 生成具有不同概率的多项式随机变量。

library(dplyr)
library(tidyr)
data %>%
  select(IDhh, Year, Educ_W, Educ_H, Educ_W1 : Educ_W5) %>% 
  nest(-IDhh) %>% 
  mutate(
    wanted_W = map(data, ~ rMultinom(t(c(.x$Educ_W1, .x$Educ_W2, .x$Educ_W3, 
                                         .x$Educ_W4, .x$Educ_W5)), 1))) %>%
  unnest()

所需的输出如下所示,其中“Wanted\u W”是新变量。

| IDhh| wanted_W|Year |Educ_W | Educ_W1| Educ_W2| Educ_W3| Educ_W4| Educ_W5| |-------:|--------:|:----|:------|---------:|---------:|---------:|---------:|---------:| | 18806| 3|1975 |3 | 0.1851884| 0.1577067| 0.4749609| 0.1394014| 0.0427427| | 2442099| 4|2010 |1 | 0.4436620| 0.0987973| 0.3296288| 0.1013606| 0.0265513| | 1351429| 3|1995 |3 | 0.0708855| 0.1023657| 0.5904598| 0.1784980| 0.0577910| | 250232| 3|1980 |5 | 0.0337913| 0.0347975| 0.2156134| 0.2315768| 0.4842209| | 1802868| 3|2005 |3 | 0.0371280| 0.0772428| 0.6054841| 0.2024385| 0.0777067| | 715077| 2|1985 |3 | 0.1149756| 0.1412112| 0.5458910| 0.1413975| 0.0565248|

3 回复 | 直到 5 年前

Mikko Marttila 7 年前

而不是打电话 Hmisc::rMultinom 一百万次(数据中的每一行一次),您可以将概率参数作为矩阵提供给函数。然后,矩阵中的每一行将定义一个不同的多项式分布。

reprex::reprex_info()
#> Created by the reprex package v0.1.1.9000 on 2018-02-09

library(dplyr)
set.seed(1)

# Generate category probabilities
n <- 1e6
unifs <- replicate(5, runif(n))
probs <- sweep(unifs, 1, apply(unifs, 1, sum), "/")
colnames(probs) <- paste0("p", seq_len(ncol(probs)))

df <- as_tibble(probs)

system.time({
  probs <- as.matrix(df %>% select(p1:p5))
  res <- df %>% 
    mutate(rcat = Hmisc::rMultinom(probs, 1))
})
#>    user  system elapsed 
#>    9.25    0.15    9.50

res
#> # A tibble: 1,000,000 x 6
#>            p1         p2         p3         p4         p5  rcat
#>         <dbl>      <dbl>      <dbl>      <dbl>      <dbl> <chr>
#>  1 0.14607852 0.07709049 0.33798110 0.22154639 0.21730349    p4
#>  2 0.12813691 0.23952958 0.11025717 0.31642808 0.20564827    p4
#>  3 0.19137423 0.24349984 0.06855848 0.23421041 0.26235703    p3
#>  4 0.30227095 0.03050219 0.27667295 0.28389810 0.10665580    p3
#>  5 0.10096040 0.03334545 0.07350112 0.38768513 0.40450791    p4
#>  6 0.32430441 0.22123172 0.13317669 0.08001760 0.24126959    p2
#>  7 0.32710720 0.14134942 0.25371663 0.20344497 0.07438178    p1
#>  8 0.21841291 0.23480314 0.25563400 0.06838794 0.22276200    p3
#>  9 0.21164692 0.19809418 0.15415735 0.15095640 0.28514514    p1
#> 10 0.02220492 0.23105648 0.35661756 0.08688459 0.30323645    p3
#> # ... with 999,990 more rows

Florian 7 年前

你可以不用 rMultinom 函数,通过生成随机均匀变量并检查其所在的区间,如所述 here :

set.seed(1234)
n=1000000
library(data.table)

# Sample data -----------------------------------------------------------
create_probs <- function(x)
{
  y = sample(1:10,x)
  y = as.list(y/sum(y))
  return(y)
}
p_dt = data.table(id=1:n)
p_dt =p_dt[,c("Educ_w1","Educ_w2","Educ_w3","Educ_w4","Educ_w5"):=create_probs(5),by=1:nrow(p_dt)]


# Function --------------------------------------------------------------
p_dt[,U:=runif(1,0,1),1:nrow(p_dt)]
p_dt = p_dt[,Educ_w:=min(which(cumsum(unlist(.SD))>U)),1:nrow(p_dt),
            .SDcols=c("Educ_w1","Educ_w2","Educ_w3","Educ_w4","Educ_w5")]

head(p_dt)

样本输出:

   id    Educ_w1    Educ_w2   Educ_w3    Educ_w4    Educ_w5          U Educ_w
1:  1 0.06666667 0.20000000 0.1666667 0.26666667 0.30000000 0.49320836      4
2:  2 0.36842105 0.05263158 0.1052632 0.26315789 0.21052632 0.54415445      4
3:  3 0.25925926 0.18518519 0.1111111 0.37037037 0.07407407 0.65840751      4
4:  4 0.29032258 0.09677419 0.3225806 0.06451613 0.22580645 0.26604797      1
5:  5 0.22222222 0.16666667 0.1111111 0.05555556 0.44444444 0.05887458      1
6:  6 0.31034483 0.17241379 0.2758621 0.20689655 0.03448276 0.98659704      5

这个功能部分在我的电脑上运行大约需要8秒钟。希望这有帮助!

Jaccar 6 年前

我在这里寻找同样问题的解决方案。老实说,我没有找到,但也许我可以给你一个更好的解决方案。

z = mapply(rmultinom, n = 1, size = 1, prob = split(probs, c(col(probs))))

所有函数都是R的内置函数, probs 是按列的,即 问题 确定一次绘制多项式的概率。结果是一个矩阵。每列都有结果(结果类别为1,其他类别为0),速度快60%