Home > voicebox > randiscr.m

randiscr

PURPOSE ^

RANDISCR Generate discrete random numbers with specified probabiities [X]=(P,N,A)

SYNOPSIS ^

function x=randiscr(p,n,a)

DESCRIPTION ^

RANDISCR Generate discrete random numbers with specified probabiities [X]=(P,N,A)

 Usage: (1) randiscr([],10)        % generate 10 uniform random binary values
        (2) randiscr(2:6,10)       % generate 10 random numbers in the range 1:5
                                     with probabilities [2 3 4 5 6]/20
        (3) randiscr([],10,'abcd') % generate a string of 10 random
                                     characters equiprobable from 'abcd'
        (4) randiscr([],-3,'abcd') % generate a string of 3 distinct random
                                     characters equiprobable from 'abcd'

 Inputs: P  vector or n-D array of probabilities (not necessarily normalized) [default = uniform]
         N  number of random values to generate: +ve=with and -ve=without replacement [default = 1]
         A  output alphabet [default = 1:length(p) or 0:1 if p is empty]

 Outputs: X  vector of not necessarily distinct values taken from alphabet A
             If P is not a vector, each row of X will contain the coordinates
             of an element of P

 The vector P is internally normalized by dividing by its sum.
 If P is an M-dimensional matrix (and A is unspecified), then X will
 have dimensions (N,M) with the corresponding indices for each dimension.

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function x=randiscr(p,n,a)
0002 %RANDISCR Generate discrete random numbers with specified probabiities [X]=(P,N,A)
0003 %
0004 % Usage: (1) randiscr([],10)        % generate 10 uniform random binary values
0005 %        (2) randiscr(2:6,10)       % generate 10 random numbers in the range 1:5
0006 %                                     with probabilities [2 3 4 5 6]/20
0007 %        (3) randiscr([],10,'abcd') % generate a string of 10 random
0008 %                                     characters equiprobable from 'abcd'
0009 %        (4) randiscr([],-3,'abcd') % generate a string of 3 distinct random
0010 %                                     characters equiprobable from 'abcd'
0011 %
0012 % Inputs: P  vector or n-D array of probabilities (not necessarily normalized) [default = uniform]
0013 %         N  number of random values to generate: +ve=with and -ve=without replacement [default = 1]
0014 %         A  output alphabet [default = 1:length(p) or 0:1 if p is empty]
0015 %
0016 % Outputs: X  vector of not necessarily distinct values taken from alphabet A
0017 %             If P is not a vector, each row of X will contain the coordinates
0018 %             of an element of P
0019 %
0020 % The vector P is internally normalized by dividing by its sum.
0021 % If P is an M-dimensional matrix (and A is unspecified), then X will
0022 % have dimensions (N,M) with the corresponding indices for each dimension.
0023 
0024 % Somewhat similar in function to RANDSRC in the comms toolbox
0025 
0026 %   Copyright (c) 2005-2012 Mike Brookes,  mike.brookes@ic.ac.uk
0027 %      Version: $Id: randiscr.m 10648 2018-07-13 15:57:20Z dmb $
0028 %
0029 %   VOICEBOX is a MATLAB toolbox for speech processing.
0030 %   Home page: http://www.ee.ic.ac.uk/hp/staff/dmb/voicebox/voicebox.html
0031 %
0032 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0033 %   This program is free software; you can redistribute it and/or modify
0034 %   it under the terms of the GNU General Public License as published by
0035 %   the Free Software Foundation; either version 2 of the License, or
0036 %   (at your option) any later version.
0037 %
0038 %   This program is distributed in the hope that it will be useful,
0039 %   but WITHOUT ANY WARRANTY; without even the implied warranty of
0040 %   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
0041 %   GNU General Public License for more details.
0042 %
0043 %   You can obtain a copy of the GNU General Public License from
0044 %   http://www.gnu.org/copyleft/gpl.html or by writing to
0045 %   Free Software Foundation, Inc.,675 Mass Ave, Cambridge, MA 02139, USA.
0046 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0047 gota=nargin>2;
0048 if nargin<1 || ~numel(p)
0049     if gota
0050         p=ones(1,length(a));
0051     else
0052         p=ones(1,2);
0053         a=(0:1)';
0054         gota=1;
0055     end
0056 end
0057 if nargin<2 || ~numel(n)
0058     n=1;
0059 end
0060 q=p(:);
0061 d=length(q);                    % size of output alphabet
0062 if n>1                          % sample with replacement
0063     dn=d+n-1;
0064     z=zeros(dn,1);               % array to hold random numbers
0065     z(1:d)=cumsum(q/sum(q));        % last value is actually overwritten in the next line
0066     z(d:end)=rand(n,1);
0067     [y,iy]=sort(z);
0068     y(iy)=(1:dn)';               % create inverse index
0069     m=zeros(dn,1);
0070     m(y(1:d-1))=1;                  % set original initial d-1 values to 1
0071     m(1)=m(1)+1;
0072     m=cumsum(m);
0073     x=m(y(d:dn));               % find the positions of the random numbers
0074 else                            % sample without replacement
0075     n=abs(n);
0076     f=(1:d)'; % list of possible outputs
0077     if n>d
0078         error('Number of output samples exceeds alphabet size');
0079     end
0080     
0081     nn=n; % number of samples remaining
0082     x=zeros(nn,1); % space for the output samples
0083     while nn>0
0084         dd=numel(f); % alphabet size
0085         qq=q(f); % alphabet probabilities
0086         if dd==1 % singleton alphabet
0087             x(n)=f;
0088         else
0089             dn=dd+nn-1;
0090             z=zeros(dn,1);               % array to hold random numbers
0091             z(1:dd)=cumsum(qq/sum(qq));        % last value is actually overwritten in the next line
0092             z(dd:dn)=rand(nn,1);
0093             [y,iy]=sort(z);
0094             y(iy)=(1:dn)';               % create inverse index
0095             m=zeros(dn,1);
0096             m(y(1:dd-1))=1;                  % set original initial d-1 values to 1
0097             m(1)=m(1)+1;
0098             m=cumsum(m);
0099             z=m(y(dd:dn));               % find the positions of the random numbers
0100             [y,iy]=sort(z);
0101             z(iy(1+find(y(1:nn-1)==y(2:nn))))=[];          % remove non-unique values
0102             k=numel(z); % number of new values
0103             x(n-nn+1:n-nn+k)=f(z);
0104             nn=nn-k;  % number of remaining samples
0105             f(z)=[]; % remove from alphabet
0106         end
0107     end
0108 end
0109 if gota
0110     x=a(x);                     % select from output alphabet
0111 elseif length(q)>length(p)      % need multiple dimensions
0112     v=x-1;
0113     s=cumprod(size(p));
0114     m=length(s);
0115     s(2:end)=s(1:end-1);
0116     s(1)=1;
0117     x=zeros(n,m);
0118     for i=m:-1:1
0119         x(:,i)=1+floor(v/s(i)); % find indices starting with the last
0120         v=rem(v,s(i));
0121     end
0122 end

Generated on Mon 06-Aug-2018 14:48:32 by m2html © 2003