File Exchange

## randp

version 3.0.01 (2.27 KB) by Jos (10584)

### Jos (10584) (view profile)

pseudorandom integers from a specified discrete distribution

Updated 02 Mar 2019

RANDP - pseudorandom integers from a specified discrete distribution

R = randp(P, N) returns an N-by-N matrix containing pseudorandom
integers drawn from a specified discrete distribution on 1:numel(P).
The distribution is specified by the relative values of P so that a
value K is present approximately "P(K)/sum(P) times in the matrix R.
All values of P should => 0, NaNs are set to 0.

The other arguments to randp specify the size of R in the same way as
matlab's own functions do: randp(P, N) returns an N-by-N matrix,
randp(P,M,N) and randp(P, [M N]) return M-by-N arrays, etc.

Examples:
% random values from [1 2 4] and a bias for 2
R = randp([1 2 0 1], 1, 100) ; % 100 values
histc(R, 1:4) % -> ~25 ~50 0 ~25

% create a random, but biased DNA sequence
C ='AGCT', P = [4 1 1 2]
DNA = C(randp(P, 1, 50))

### Cite As

Jos (10584) (2020). randp (https://www.mathworks.com/matlabcentral/fileexchange/8891-randp), MATLAB Central File Exchange. Retrieved .

Jeff Miller

Ali Komai

### Ali Komai (view profile)

Marcelo Fernandes

### Marcelo Fernandes (view profile)

Thanks for this file! Works great!

Martin

### Martin (view profile)

What an elegant piece of code! Thanks for this.

Liber Eleutherios

### Liber Eleutherios (view profile)

I've tried this small experiment:

N = 100000; % sample size
w = 1:100; % quite steep distribution
w = w / sum(w); % distribution (pdf)
Y_std = sqrt(w .* (1 - w)); % Standard deviations associated to each probability

%%%%%%%%%%%%%%%%%

Y1 = randp(w, N, 1); % Try Jos function
Y1_w = histc(Y1, [1:100]); % absolute frequencies table
Y1_w = Y1_w / N; % transform into relative frequencies
Y1_w_std = sqrt(Y1_w .* (1 - Y1_w)); % Empirical standard deviations associated to each probability

%%%%%%%%%%%%%%%%%

Y2 = gDiscrPdfRnd(w, N, 1); % Try Gianluca Dorini's function
Y2_w = histc(Y2, [1:100]); % absolute frequencies table
Y2_w = Y2_w / N; % transform into relative frequencies
Y2_w_std = sqrt(Y2_w .* (1 - Y2_w)); % Empirical standard deviations associated to each probability

%%%%%%%%%%%%%%%%%

% Check if empirical standard deviations are coherent with theretical ones
plot(Y_std - Y1_w_std', '.b') % plot empirical differences for Jos function
hold on
plot(Y_std - Y2_w_std', '.r') % plot empirical differences for Gianluca Dorini's function

%%%%%%%%%%%%%%%%%

It looks as though both functions statistically behave very well. They are both excellent functions. Thank you.

Liber Eleutherios

### Liber Eleutherios (view profile)

Wow, it really is a nice function, very fast and efficient. It is exactly what I was looking for. Thank you.

Carlos Baiz

### Carlos Baiz (view profile)

Works great!

Very useful when dealing with user-defined PDFs.

urs (us) schwarz

nicely coded snippet
1) why not put the validity test for P BEFORE the (potentially consuming) generation of all the random numbers...
2) why not add a pointer to RANDSAMPLE for those who own the stats tbx...
us