Hi,


I have a set of N objects randomly distributed in a 2D physical space. Each object (i) generates a bernoulli random number (0 or 1) based on a marginal probability Pr(xi = 1) = p. These objects a correlated by physical distance. The closer the objects are, the larger their correlation is.


E.g. If objects i and j are co-located, they are expected to generate correlated results. For Example, if P(Xi=1)= 0.6 and P(Xj=1)=0.3 they would produce something like:


Xi= 0 1 0 1 1 1 0 1 0 1


Xj= 0 1 0 0 0 1 0 1 0 0


Such that Pr(Xi|Xj)=1


On the other hand if i and j are distant they would produce uncorrelated results such that
Pr(Xi|Xj)=Pr(Xi)


I have tried to use some of the packages in Matlab (Sampling from multivariate correlated binary and poisson random variables) and R (bindata) but I could not produce an acceptable correlation matrix.


Any ideas how I can produce an acceptable correlation matrix?


BTW, I have checked the following earlier posts [discrete stochastic process: exponentially correlated Bernoulli?](https://mathoverflow.net/questions/18268/discrete-stochastic-process-exponentially-correlated-bernoulli)


and


[Constructing Bernoulli random variables with prescribed correlation](https://mathoverflow.net/questions/19406/constructing-bernoulli-random-variables-with-prescribed-correlation)


But I am not sure how I can relate to them.


Thanks



Generating Bernoulli Correlated Random Variables with Space Decaying Correlations
Here's a suggestion:


Define a non-negative decreasing function $w(r)$ measuring interaction strength. Given each 
object its own independent $N(0,1)$ random variable $N\_i$. Now set
$$
Y\_i=\frac{\sum\_{j}w(\|x\_i-x\_j\|)N\_j}{\sqrt{\sum\_j w(\|x\_i-x\_j\|)^2}},
$$
where $x\_i$ denotes the location of the $i$th object.


Then the $Y\_i$ are correlated $N(0,1)$ random variables. If two objects are co-located the normal random variables agree. 


Finally set $t\_i=\Phi^{-1}(p\_i)$ (i.e. $\mathbb P(N < t\_i)=p\_i$) and set $X\_i=1$ if $Y\_i < p\_i$ and 0 otherwise. 


With this setup you can write down the covariance of $Y\_i$ and $Y\_k$ explicitly: it's just $$
\text{Cov}(Y\_i,Y\_k)=\frac{\sum\_j w(\|x\_i-x\_j\|)w(\|x\_k-x\_j\|)}
{\sqrt{\sum\_j w(\|x\_i-x\_j\|)^2\sum\_j w(\|x\_k-x\_j\|)^2}}.
$$


If you write this as $\cos\theta\_{ik}$ then you can write the covariance of $X\_i$ and $X\_k$ as an integral: 
$$
1/(2\pi)\int\_{
x < t\_1\;,\;
cos\theta\_{ik}x+\sin\theta\_{ik}y < t\_2}
e^{-(x^2+y^2)/2}\,dxdy-p\_ip\_k.
$$


