时间戳

首页 > 解决方案 > 模拟后两个向量的比较

r - 模拟后两个向量的比较

问题描述

我想应用拒绝采样方法来模拟Y=(Y_1, Y_2)来自单位圆盘D = { (X_1 , X_2) \in R^2: \sqrt{x^2_1 + x^2_2} ≤ 1}的X = (X_1 , X_ 2)均匀分布的随机向量，即正方形中均匀分布的随机向量S = [−1, 1]^2和联合密度f(y_1,y_2) = \frac{1}{\pi} 1_{D(y_1,y_2)}.

在拒绝方法中，如果 ，我们通常接受一个样本f(x) \leq C * g(x)。我正在使用以下代码：

x=runif(100,-1,1)
y=runif(100,-1,1)

d=data.frame(x=x,y=y)
disc_sample=d[(d$x^2+d$y^2)<1,]
plot(disc_sample)

我有两个问题：

{使用上面的代码，逻辑上，大小d应该大于大小，disc_sample但是当我调用它们时，我看到它们每个都有 100 个元素。这怎么可能。为什么尺寸相同。}这部分已解决，感谢下面的评论。

现在的问题

此外，我如何重新编写我的代码，以提供使 100 个样本符合条件所需的样本总数。即在我得到所需的 100 个样本之前给我拒绝的样本数量？

感谢r2evans的回答，但我希望写一些更简单的东西，一个 while 循环将所有可能的样本存储在矩阵或数据帧而不是列表中，然后从该数据帧调用只是样本遵循条件。我修改了答案中的代码，没有使用列表，也没有 sapply 函数，但它没有给出所需的结果，它只产生一行。

i=0
samps <- data.frame()
goods <- data.frame()
nr <- 0L
sampsize <- 100L
needs <- 100L
while (i < needs) {
  samps <- data.frame(x = runif(1, -1, 1), y = runif(1, -1, 1))
  goods <- samps[(samps$x^2+samps$y^2)<1, ]
i = i+1
}

我也想过这个：

i=0
j=0
samps <- matrix()
goods <- matrix()
needs <- 100

while (j < needs) {
  samps[i,1] <- runif(1, -1, 1)
  samps[i,2] <- runif(1, -1, 1)
  if (( (samps[i,1])**2+(samps[i,2])**2)<1){
  goods[j,1] <- samps[i,1]
  goods[j,2] <- samps[i,2]
}
else{
i = i+1
}
}

但它不工作。

对于修改代码的任何帮助，我将不胜感激。

标签： rsimulationprobabilitysampling

解决方案

至于你的第二个问题......你不能重新编写你的代码来准确地知道获得（至少）100 个结果组合需要多少。您可以使用while循环并连接结果，直到您至少有 100 个这样的行，然后截断超过 100 个的行。因为分段（按比例）使用熵是“昂贵的”，您可能更喜欢总是高估您需要的行和一次抓住所有。

（编辑以减少基于作业限制的“复杂性”。）

set.seed(42)
samps <- vector(mode = "list")
goods <- vector(mode = "list")
nr <- 0L
iter <- 0L
sampsize <- 100L
needs <- 100L
while (nr < needs && iter < 50) {
  iter <- iter + 1L
  samps[[iter]] <- data.frame(x = runif(sampsize, -1, 1), y = runif(sampsize, -1, 1))
  rows <- (samps[[iter]]$x^2 + samps[[iter]]$y^2) < 1
  goods[[iter]] <- samps[[iter]][rows, ]
  nr <- nr + sum(rows)
}
iter                # number of times we looped
# [1] 2
out <- head(do.call(rbind, goods), n = 100)
NROW(out)
# [1] 100
head(out) ; tail(out)
#            x          y
# 1  0.8296121  0.2524907
# 3 -0.4277209 -0.5668654
# 4  0.6608953 -0.2221099
# 5  0.2834910  0.8849114
# 6  0.0381919  0.9252160
# 7  0.4731766  0.4797106
#               x          y
# 221 -0.65673577 -0.2124462
# 231  0.08606199 -0.7161822
# 251 -0.37263236  0.1296444
# 271 -0.38589120 -0.2831997
# 28  -0.62909284  0.6840144
# 301 -0.50865171  0.5014720

推荐阅读