いくつかの二項分布から相関乱数を生成する方法を見つけようとしています。
( を使用して) 正規分布でそれを行う方法は知ってMASS::mvrnorm
いますが、二項応答に適用できる関数が見つかりませんでした。
いくつかの二項分布から相関乱数を生成する方法を見つけようとしています。
( を使用して) 正規分布でそれを行う方法は知ってMASS::mvrnorm
いますが、二項応答に適用できる関数が見つかりませんでした。
copula
パッケージを使用して相関ユニフォームを生成し、qbinom
関数を使用してそれらを二項変数に変換できます。簡単な例を次に示します。
library(copula)
tmp <- normalCopula( 0.75, dim=2 )
x <- rcopula(tmp, 1000)
x2 <- cbind( qbinom(x[,1], 10, 0.5), qbinom(x[,2], 15, 0.7) )
x2
これで、相関する2つの二項変数を表す2つの列を持つ行列になります。
n 回の試行と各試行での成功確率 p を持つ二項変数は、n 回のベルヌーイ試行の合計と見なすことができ、それぞれが成功確率 p を持ちます。
同様に、目的の相関 r を持つベルヌーイ変量のペアを合計することで、相関のある二項変量のペアを作成できます。
require(bindata)
# Parameters of joint distribution
size <- 20
p1 <- 0.5
p2 <- 0.3
rho<- 0.2
# Create one pair of correlated binomial values
trials <- rmvbin(size, c(p1,p2), bincorr=(1-rho)*diag(2)+rho)
colSums(trials)
# A function to create n correlated pairs
rmvBinomial <- function(n, size, p1, p2, rho) {
X <- replicate(n, {
colSums(rmvbin(size, c(p1,p2), bincorr=(1-rho)*diag(2)+rho))
})
t(X)
}
# Try it out, creating 1000 pairs
X <- rmvBinomial(1000, size=size, p1=p1, p2=p2, rho=rho)
# cor(X[,1], X[,2])
# [1] 0.1935928 # (In ~8 trials, sample correlations ranged between 0.15 & 0.25)
望ましい相関係数を共有するさまざまな結合分布が多数あることに注意することが重要です。のシミュレーション メソッドrmvBinomial()
はそれらの 1 つを生成しますが、それが適切かどうかは、データを生成するプロセスによって異なります。
同様の質問に対するこのRヘルプの回答に記載されているように(その後、アイデアをより詳細に説明します):
2 変量正規 (与えられた平均と分散) は相関係数によって一意に定義されますが、これは 2 変量二項には当てはまりません。