给定二元离散分布的随机样本

Question

假设我有一个二元离散分布，即i=1，...n和j=1，...m的概率值P(X=i，Y=j)表。如何从这样的分布中生成随机样本(X_k，Y_k)，k=1，...N？也许有一个现成的R函数，比如：

sample(100,prob=biprob)

其中biprob是2维矩阵？

下面是一种直观的采样方法。假设我们有一个data.frame

dt=data.frame(X=x,Y=y,P=pij)

X和y从哪里来

expand.grid(x=1:n,y=1:m)

和pij是P(X=i，Y=j)。

然后，我们通过以下方式获得大小为N的样本(Xs，Ys)：

set.seed(1000) 
Xs <- sample(dt$X,size=N,prob=dt$P)
set.seed(1000)
Ys <- sample(dt$Y,size=N,prob=dt$P)

我使用set.seed()来模拟“双变量”。直觉上，我应该得到与我需要的东西相似的东西。不过，我不确定这是不是正确的方式。因此出现了这样的问题:)

另一种方法是使用Gibbs抽样，边缘分布很容易计算。

我试着用谷歌搜索，但没有找到真正相关的内容。

Answer 1

你就快到了。假设您有包含x、y和pij值的数据框dt，只需对这些行进行采样！

dt <- expand.grid(X=1:3, Y=1:2)
dt$p <- runif(6)
dt$p <- dt$p / sum(dt$p)  # get fake probabilities
idx <- sample(1:nrow(dt), size=8, replace=TRUE, prob=dt$p)
sampled.x <- dt$X[idx]
sampled.y <- dt$Y[idx]

Answer 2

我不明白为什么你要关心它是二元的。概率和为1，并且结果是离散的，因此您只是从categorical distribution中采样。唯一的区别是，您使用行和列而不是单个位置对观察结果进行索引。这只是一个符号。

因此，在R中，您可以通过重塑数据并从分类分布中采样来轻松地从您的分布中采样。可以使用rmultinom和which选择索引，或者像Aniko建议的那样，使用sample对整形后的数据行进行采样。一些记账员可以处理你的具体情况。

这里有一个解决方案：

library(reshape)

# Reshape data to long format.
data <- matrix(data = c(.25,.5,.1,.4), nrow=2, ncol=2)
pmatrix <- melt(data)

# Sample categorical n times.
rcat <- function(n, pmatrix) {
    rows <- which(rmultinom(n,1,pmatrix$value)==1, arr.ind=TRUE)[,'row']
    indices <- pmatrix[rows, c('X1','X2')]
    colnames(indices) <- c('i','j')
    rownames(indices) <- seq(1,nrow(indices))
    return(indices)
}

rcat(3,pmatrix)

这将从矩阵中返回3个随机抽样值，报告行和列的i和j：

  i j
1 1 1
2 2 2
3 2 2