2 Probability Theory 2.2 Exact Formula for Asymptotic Convergence of Fourier Transform of Uniform Random Variables 2.4 Distributions of Functions of Random Variables and Matrices

2.3 Generating and Sampling Two Signed Bernoulli Random Variables with Given Correlation

The following first constructs two signed $\{-1,+1\}$ Bernoulli variables $Y_{1},Y_{2}$ such that $p(Y_{1}=-1)=p_{1},p(Y_{2}=-1)=p2,$ and $E[Y_{1}Y_{2}]=r.$ Constraints: $0<p1,p2<1,-1<r<1.$

Next, it draws $N$ realizations of the two variables in $Y(m,:),m=1,2.$ The sample pdf and cross correlation converge as $N\to\infty.$

2.3.1 MATLAB Reference Code

To generate $N$ realizations of variables with mean $0,$ variance $1,$ and correlation $0.3$ using the code below, use: Y = signed_bernoulli2(N,0.5,0.5,0.3).

⬇

function Y = signed_bernoulli2(N,p1,p2,r)

% Y = signed_bernoulli2(N,p1,p2,r)

% The function first constructs two(2) signed {-1, +1} Bernoulli variables

% Y_1, Y_2 such that p(Y_1 = -1) = p1, p(Y_2 = -1) = p2, and E[Y_1Y_2] = r.

% Constraints: 0 < p1,p2 < 1, -1 < r < 1.

% Next, it draws N realizations of the two variables in Y(m,:), m = 1,2.

% The sample pdf and cross correlation converge as N -> infinity (obviously).

% Example: generate N realizations of variables with mean 0, variance 1, and

% correlation 0.3: Y = signed_bernoulli2(N,0.5,0.5,0.3).

% Author: aravindh.krishnamoorthy@fau.de, on: 26.10.2017.

% Basic checks and initialization

assert (0 < p1) ;

assert (p1 < 1) ;

assert (0 < p2) ;

assert (p2 < 1) ;

assert (-1 < r) ;

assert (r < 1) ;

Y = zeros(2,N) ;

% Resolve the joint pdf given by

% / p(1) if (+1, +1)

% p(Y_1, Y_2) = | p(2) if (+1, -1)

% | p(3) if (-1, +1)

% \ p(4) if (-1, -1)

% Conditions (i.e. how matrices A and b are generated):

% row 1: p(1) + p(2) + p(3) + p(4) = 1

% row 2: p(3) + p(4) = p1 -> marginal for Y_1

% row 3: p(2) + p(4) = p2 -> marginal for Y_2

% row 4: p(1) - p(2) - p(3) + p(4) = r -> cross correlation E[Y_1 Y_2]

A = [1,1,1,1;0,0,1,1;0,1,0,1;1,-1,-1,1] ;

b = [1;p1;p2;r] ;

p = pinv(A)*b ;

assert(all(p >= 0), sprintf('The given p1=%f and p2=%f cannot result in a correlation r=%f with the given algorithm.\np1=p2=0.5 is a safe value.',p1,p2,r)) ;

% Now that we have the joint pdf in "p" we can draw from this distribution

% as follows (based on conditional pdfs and inverse transform sampling):

% 1. Draw samples of Y_1 from Signed_Bernoulli(p1) in a standard fashion using

% an auxiliary uniformly distributed variable U_1 = U(1,:).

% 2. If the n-th sample of Y_1 is "-1", then draw for Y_2 from the conditional pdf

% p(Y_2|Y_1 = -1) using a second axiliary uniformly distributed

% variable U_2 = U(2,:).

% 2. If the n-th sample of Y_1 is "+1", then draw for Y_2 from the conditional pdf

% p(Y_2|Y_1 = +1) using the second axiliary uniformly distributed

% variable U_2 = U(2,:).

% The conditional pdfs are given as follows:

% p(Y_2|Y_1 = -1) = p(4)/p1 for -1, and 1-(p(4)/p1) for +1.

% p(Y_2|Y_1 = +1) = p(2)/(1-p1) for -1, and 1-(p(2)/(1-p1)) for +1.

% Note: for the rather straightforward proofs for the above steps, please contact

% the author.

% Auxiliary variable U

U = rand(2,N) ;

% Draw samples from p(Y_1)

Y(1,:) = (-1).^(U(1,:) < p1) ;

% Indices to help choose the conditional pdf

I1 = Y(1,:) == -1 ;

I2 = Y(1,:) == +1 ;

% Draw samples from the conditional pdfs

Y(2,I1) = (-1).^(U(2,I1) < p(4)/p1) ;

Y(2,I2) = (-1).^(U(2,I2) < p(2)/(1-p1)) ;

The full program can also be downloaded from signed_bernoulli2(N,p1,p2,r) - File Exchange - MATLAB Central (mathworks.com).

2.3.2 Version History

1.

First published: 23rd Nov. 2017 on aravindhk-math.blogspot.com
2.

Modified: 17th Dec. 2023 – Style updates for LaTeX
3.

Modified: 30th Dec. 2023 – Minor updates to text