如何生成具有共同相關性的伯努利隨機變量ρρrho?

假設我想生成伯努利隨機變量，例如：

對所有人.

我想知道我可以使用什麼方法（我將嘗試在 R 中執行此操作）。我正在考慮使用 Couplas，但有人可以指導我如何實現嗎？謝謝。

因為這個相關矩陣是如此對稱，我們可能會嘗試用對稱分佈來解決這個問題。

以下是在改變相關性方面提供足夠靈活性的最簡單方法之一。給定變量，定義集合上的分佈維二進制向量通過分配概率到, 概率到，並分配剩餘概率同樣在正好有一個的向量; 因此，這些中的每一個都有概率. 請注意，這一系列分佈僅取決於一個參數.

**從這些分佈之一很容易模擬：**輸出一個具有概率的零向量，輸出一個有概率的向量，否則從單位矩陣。

的所有組件是同分佈的伯努利變量。 它們都有共同的參數

計算協方差和通過觀察它們都可以相等只有當所有組件都, 從何而來

這將互相關確定為

給定和（這是任何可能的相關性的範圍-variate 隨機變量），存在唯一解之間和.

模擬證明了這一點。 從一組開始等距的值，對應的值被計算（對於這種情況) 並用於生成的獨立值. 這計算相關係數並將其繪製在垂直軸上。協議很好。

我對以下值進行了一系列此類模擬之間和, 效果比較好。

這種方法的一般化（即，允許兩個、三個或 …同時等於) 將提供更大的靈活性，在此解決方案中由下式確定. 這將這里相關的想法與完全一般的想法結合在一起https://stats.stackexchange.com/a/285008/919中描述的解決方案。

---

以下R代碼具有p計算功能從和並在其主循環中展示了相當有效的模擬機制。

#
# Determine p(All zeros) = p(All ones) from rho and d.
#
p <- function(rho, d) {
 if (rho==1) return(1/2)
 if (rho <= -1/(d-1)) return(0)
 if (d==2) return((1+rho)/4)
 b <- d-2
 (4 + 2*b + b^2*(1-rho) - (b+2)*sqrt(4 + b^2 * (1-rho)^2)) / (2 * b^2 * (1-rho))
}
#
# Simulate a range of correlations `rho`.
#
d <- 8       # The number of variables.
n.sim <- 1e4 # The number of draws of X in the simulation.

rho.limits <- c(-1/(d-1), 1)
rho <- seq(rho.limits[1], rho.limits[2], length.out=21)
rho.hat <- sapply(rho, function(rho) {
 #
 # Compute the probabilities from rho.
 #
 qd <- q0 <- p(rho, d)
 q1 <- (1 - q0 - qd)
 #
 # First randomly select three kinds of events: all zero, one 1, all ones.
 #
 u <- sample.int(3, n.sim, prob=c(q0,q1,qd), replace=TRUE)
 #
 # Conditionally, when there is to be one 1, uniformly select which 
 # component will equal 1.
 #
 k <- diag(d)[, sample.int(d, n.sim, replace=TRUE)]
 #
 # When there are to be all zeros or all ones, make it so.
 #
 k[, u==1] <- 0
 k[, u==3] <- 1
 #
 # The simulated values of X are the columns of `k`. Return all d*(d-1)/2 correlations.
 #
 cor(t(k))[lower.tri(diag(d))]
})
#
# Display the simulation results.
#
plot(rho, rho, type="n",
    xlab="Intended Correlation",
    ylab="Sample Correlation", 
    xlim=rho.limits, ylim=rho.limits,
    main=paste(d, "Variables,", n.sim, "Iterations"))

abline(0, 1, col="Red", lwd=2)

invisible(apply(rho.hat, 1, function(y) 
 points(rho, y, pch=21, col="#00000010", bg="#00000004")))

引用自：https://stats.stackexchange.com/questions/318759