Genera numeri casuali correlati dalle distribuzioni binomiali in R

Sto cercando di trovare un modo per generare numeri casuali correlati da diverse distribuzioni binomiali.

So come farlo con le normali distribuzioni (usando mvrnorm), ma non ho trovato una funzione applicabile a quelle binomiali.

È ansible generare uniformi correlate utilizzando il pacchetto copula , quindi utilizzare la funzione qbinom per convertirle in variabili binomiali. Ecco un rapido esempio:

 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) ) 

Ora x2 è una matrice con le 2 colonne che rappresentano 2 variabili binomiali che sono correlate.

Una variabile binomiale con n prove e probabilità p di successo in ogni prova può essere vista come la sum di n prove di Bernoulli ciascuna avendo anche probabilità di successo.

Allo stesso modo, è ansible build coppie di variabili binomiali correlate sumndo le coppie di variabili di Bernoulli aventi la correlazione desiderata 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) 

È importante notare che ci sono molte diverse distribuzioni congiunte che condividono il coefficiente di correlazione desiderato. Il metodo di simulazione in rmvBinomial() produce uno, ma indipendentemente dal fatto che sia appropriato dipenderà dal processo che genera i dati.

Come notato in questa risposta R-help a una domanda simile (che poi continua a spiegare l'idea in modo più dettagliato):

mentre una normale bivariata (dati mezzi e varianze) è definita in modo univoco dal coefficiente di correlazione, questo non è il caso per un binomio bivariato