Take the 2-minute tour ×
MathOverflow is a question and answer site for professional mathematicians. It's 100% free, no registration required.

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?

and

Constructing Bernoulli random variables with prescribed correlation

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

Thanks

share|improve this question

1 Answer 1

up vote 1 down vote accepted

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. $$

share|improve this answer
    
If i understood you correctly, what you are saying is to simply use the Cov(Yi,Yk) above as the covariance matrix and plug it in the bernoulli generator. If so, then Unfortunately, I tried that but I keep getting an Unacceptable correlation matrix, which seams to be not positive definite. –  alandalusi Apr 26 '11 at 21:44
    
No you didn't understand me correctly. Here is my concrete suggestion. (1) Compute $t_i=\Phi^{-1}(p_i)$; (2) Compute the matrix $Cov(Y_i,Y_k)$ as above (3) Use a multivariate normal generator to build some $Y_i$'s. (4) Set $X_i=1$ if $Y_i<t_i$ and 0 otherwise. You could compute the covariance of the Bernoulli's (using the nasty integral formula I wrote down), but if your purpose is to just generate the random #s, then the procedure I described will work just as well. –  Anthony Quas Apr 26 '11 at 22:59
    
Acually, It is positive definite, but does not have a dichotomized Gaussian distribution for the correlation matrix. –  alandalusi Apr 26 '11 at 23:03
    
OK, I think I understand your suggestion now. I will test it right away. –  alandalusi Apr 26 '11 at 23:32
    
It works!, but I noticed that the Correlation (or Conditional Probabilities) $Pr(X{_i}=1|X{_k}=1)$ between any fixed pair (i and k) slightly change as we add more objects into the physical space. So, I changed $Cov(Y{_i},Y{_k})$ to be $w(||x{_i}-x{_k}||)$ and it worked with fixed conditional probabilities, even if I add more objects. I am still testing it, but does this make sense to you? For your information, my distance functions is $e^{(-\frac{eucaliandistance}{dcorr})}$ –  alandalusi Apr 27 '11 at 16:33

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.