在R中设置随机采样时间

disho6za  于 2023-04-18  发布在  其他
关注(0)|答案(2)|浏览(80)

我正在使用R for NONMEM创建模拟数据。
我想知道如何为34个样本创建3个样本时间的随机集。
下面是我的代码:

nsubjects <- 34

times <- rep((c(runif(1, 0.5, 1), runif(1, 12, 48), runif(1, 47, 49))), nsubjects)

现在,唯一的问题是这组三个一组的样本没有个体间的变异性,在真实的世界中,并不是每个人的每个样本都是在完全相同的时间采集的,这就是为什么我使用runif命令来设置条件,第一个样本必须在0.5-1小时之间采集,第二个样本必须在12-48小时之间采集,第三个样本必须在47-49小时之间采集。
我不在乎这些数据是否在数组中,它们只需要位于一个向量中(所以没有矩阵,我不想要三个采样时间列,因为我将数据附加到ID中。
每组3个样本必须在其预期运行的参数范围内(最小值-最大值),应有3 * 34个样本,每个样本都有其自己的预期随机值。

hfwmuf9z

hfwmuf9z1#

R的伪RNG是矢量化的,定义了下限和上限,一个对runif的调用就可以完成这项工作。
在下面的代码中,我将矢量化的结果与Ric的注解中的结果进行了比较。

nsubjects <- 34L

set.seed(2023)    # make the results reproducible
times1 <- replicate(nsubjects, c(runif(1, 0.5, 1), runif(1, 12, 48), runif(1, 47, 49)))

set.seed(2023)    # reproduce the results above
lo <- c(0.5, 12, 47)
hi <- c(1, 48, 49)
times2 <- replicate(nsubjects, runif(3, lo, hi))

all(times1 == times2)
#> [1] TRUE

identical(times1, times2)
#> [1] TRUE

创建于2023-04-06带有reprex v2.0.2
你也可以用length(lo)(或hi)替换runif调用中的常量3

times2 <- replicate(nsubjects, runif(length(lo), lo, hi))

由于您需要的是向量,而不是矩阵,假定R的矩阵是列一阶的,只需将上面的结果强制为向量c()

times2 <- c(times2)

或直接

times <- c(replicate(nsubjects, runif(3, lo, hi)))

完整的代码是

nsubjects <- 34L

set.seed(2023)    # make the results reproducible
lo <- c(0.5, 12, 47)
hi <- c(1, 48, 49)
times2 <- c(replicate(nsubjects, runif(length(lo), lo, hi)))
times2
#>   [1]  0.7333070 24.0668743 47.3256351  0.6980600 13.0941023 47.2417697
#>   [7]  0.7130828 34.2428837 47.5264165  0.7381619 43.0434840 47.2977573
#>  [13]  0.5902150 47.9738392 48.6834824  0.5713245 24.4145419 48.7792509
#>  [19]  0.6599241 37.5245700 48.2563034  0.7190863 24.2362570 48.8661852
#>  [25]  0.8564772 34.7106674 48.4928989  0.5557251 21.8383865 47.1970362
#>  [31]  0.8787647 35.8893157 48.4487222  0.7294733 36.2783489 48.4010272
#>  [37]  0.9075361 36.8825336 48.4328247  0.9452899 24.2506878 48.3371257
#>  [43]  0.6921270 20.4640812 48.7769941  0.9007258 33.4713985 48.9540355
#>  [49]  0.6640525 40.9864810 47.0646308  0.7241449 21.5303418 47.8185350
#>  [55]  0.5972984 31.5655505 48.8702074  0.8870188 47.8877970 48.6452517
#>  [61]  0.6958527 25.8689813 48.7982756  0.6594265 19.5127854 48.1283543
#>  [67]  0.7410422 14.3441268 47.5225547  0.6651678 23.9003085 47.1847066
#>  [73]  0.8469968 16.1936265 48.7411946  0.5643353 42.0736070 47.4753838
#>  [79]  0.7252093 42.0347502 48.3315939  0.7720596 35.8425095 48.7470245
#>  [85]  0.8341940 19.6876291 47.8381487  0.7313880 34.2771132 47.9674524
#>  [91]  0.8781355 19.9780284 48.9823589  0.9898648 42.9250703 47.4830284
#>  [97]  0.8156663 38.7856431 48.5488074  0.8374183 38.8840239 48.9481844

创建于2023-04-06带有reprex v2.0.2

编辑

这里有一个避免merge的解决方案。创建一个临时列count,其中包含所需的ID重复次数。uncount该列的数据。运行上面的c(replicate(etc))代码来创建新列。
末尾的结果显示指定范围内的数字。

suppressPackageStartupMessages({
  library(dplyr)
  library(tidyr)
})

nsubjects <- 34L
df1 <- data.frame(ID = seq.int(nsubjects))

lo <- c(0.5, 12, 47)
hi <- c(1, 48, 49)

set.seed(2023)    # make the results reproducible
df1 %>%
  mutate(count = length(lo)) %>%
  uncount(count) %>%
  mutate(TIME = c(replicate(nsubjects, runif(length(lo), lo, hi)))) -> df2

head(df2, n = 10)
#>    ID       TIME
#> 1   1  0.7333070
#> 2   1 24.0668743
#> 3   1 47.3256351
#> 4   2  0.6980600
#> 5   2 13.0941023
#> 6   2 47.2417697
#> 7   3  0.7130828
#> 8   3 34.2428837
#> 9   3 47.5264165
#> 10  4  0.7381619

创建于2023-04-06使用reprex v2.0.2

ih99xse1

ih99xse12#

请注意,runif是矢量化的。只需执行以下操作:

runif(34 * 3,  c(0.5, 12, 47),c(1, 48, 49))

这个向量的头是:

[1]  0.9758125 12.1717997 48.5416910  0.5152847 21.6146284 47.9904000

相关问题