5.0

5.0 | 2 ratings Rate this file 24 Downloads (last 30 days) File Size: 6.9 KB File ID: #35797
image thumbnail

Generate Random Numbers from a 2D Discrete Distribution

by

 

21 Mar 2012 (Updated )

Random numbers from any 2D discrete probability distribution, at any resolution.

| Watch this File

File Information
Description

Tristan Ursell
2D Random Number Generator for a Given Discrete Distribution
March 2012
[x0,y0]=pinky(Xin,Yin,dist_in,varargin);
'Xin' is a vector specifying the equally spaced values along the x-axis.

'Yin' is a vector specifying the equally spaced values along the y-axis.

'dist_in' (dist_in > 0) is a matrix with dimensions length(Yin) x length(Xin), whose values specify a 2D discrete probability distribution. The distribution does not need to be normalized.

'res' (res > 1) is a multiplicative increase in the resolution of chosen random number, using cubic spline interpolation of the values in 'dist_in'. Using the 'res' option can significantly slow down the code, due to the computational costs of interpolation, but allows one to generate more continuous values from the distribution.

[x0,y0] is the output doublet of random numbers consistent with dist_in.

The 'gendist' function required by this script, is included in this m-file.

Example with an anisotropic Gaussian at native resolution:

 Xin=-10:0.1:10;
 Yin=-5:0.1:5;

 Xmat = ones(length(Yin),1)*Xin;
 Ymat = Yin'*ones(1,length(Xin));

 D1=1;
 D2=3;
 Dist=exp(-Ymat.^2/(2*D1^2)-Xmat.^2/(2*D2^2));

 N=10000;

 vals=zeros(2,N);
 for i=1:N
    [vals(1,i),vals(2,i)]=pinky(Xin,Yin,Dist);
 end

 figure;
 hold on
 imagesc(Xin,Yin,Dist)
 colormap(gray)
 plot(vals(1,:),vals(2,:),'r.')
 xlabel('Xin')
 ylabel('Yin')
 axis equal tight
 box on

Example with multiple peaks at 10X resolution:

 Xin=-10:0.1:10;
 Yin=-5:0.1:5;

 Xmat = ones(length(Yin),1)*Xin;
 Ymat = Yin'*ones(1,length(Xin));

 D1=0.5;
 D2=1;
 Dist=exp(-(Ymat-3).^2/(4*D2^2)-(Xmat+6).^2/(2*D1^2))+exp(-(Ymat+2).^2/(2*D1^2)-(Xmat+1).^2/(2*D2^2))+exp(-(Ymat-1).^2/(2*D1^2)-(Xmat-2).^2/(2*D2^2));

 N=10000;
 res=10;
 vals=zeros(2,N);
 for i=1:N
    [vals(1,i),vals(2,i)]=pinky(Xin,Yin,Dist,res);
 end

 figure;
 hold on
 imagesc(Xin,Yin,Dist)
 colormap(gray)
 plot(vals(1,:),vals(2,:),'r.')
 xlabel('Xin')
 ylabel('Yin')
 axis equal tight
 box on

Required Products MATLAB
MATLAB release MATLAB 8.3 (R2014a)
MATLAB Search Path
/
Tags for This File   Please login to tag files.
Please login to add a comment or rating.
Comments and Ratings (3)
17 Oct 2014 Tristan Ursell

Good catch -- I just put a check in the code to make sure the interpolated values are always positive.

01 Oct 2014 M S

Excellent code, and amazingly, the only FEX submission for random sampling from 2D distributions.

I did run into a little bug. When using interpolation (res>1), I occasionally get "Error: All elements of first argument, P, must be positive." thrown by gendist(). This looks like it comes from interpolation giving negative numbers. The only hint at the culprit is that the joint PDFs I'm sending to your function have very very small numbers, (e.g., 1.4822e-323). When more reasonable numbers are used, the interpolation is fine.

26 Apr 2012 Filip Trönnberg

Wonderful function!!

Updates
17 Oct 2014

fixed negative values due to interpolation

17 Oct 2014

bug fix

Contact us