V_RANDVEC Generate real or complex GMM/lognormal random vectors X=(N,M,C,W,MODE)
generates a random matrix of size (|n|,p) where p is the maximum
dimension of M or C (see note below about row versus column vectors)
Inputs: N is the number of points to generate
M(K,P) is the mean vectors (one row per mixture)
C(K,P) are diagonal covariances (one row per mixture)
or C(P,P,K) are full covariance matrices (one per mixture)
W(K) are the mixture weights (or omit if all mixtures have equal weight)
MODE character string specifying options:
g = real gaussian [default]
c = complex gaussian
l = lognormal
Outputs: X(N,P) is the output data
Note this routine generates row vectors such that E((x-m)'*(x-m)) = C = cov(x). If
Alternatively x=(n,m',c,w,mode)' will generate column vectors satisfying E((x-m)*(x-m)') = C = cov(x').0001 function x=v_randvec(n,m,c,w,mode) 0002 %V_RANDVEC Generate real or complex GMM/lognormal random vectors X=(N,M,C,W,MODE) 0003 % generates a random matrix of size (|n|,p) where p is the maximum 0004 % dimension of M or C (see note below about row versus column vectors) 0005 % Inputs: N is the number of points to generate 0006 % M(K,P) is the mean vectors (one row per mixture) 0007 % C(K,P) are diagonal covariances (one row per mixture) 0008 % or C(P,P,K) are full covariance matrices (one per mixture) 0009 % W(K) are the mixture weights (or omit if all mixtures have equal weight) 0010 % MODE character string specifying options: 0011 % g = real gaussian [default] 0012 % c = complex gaussian 0013 % l = lognormal 0014 % 0015 % Outputs: X(N,P) is the output data 0016 % 0017 % Note this routine generates row vectors such that E((x-m)'*(x-m)) = C = cov(x). If 0018 % Alternatively x=(n,m',c,w,mode)' will generate column vectors satisfying E((x-m)*(x-m)') = C = cov(x'). 0019 0020 % Bugs/suggestions 0021 % (1) New mode 'x' to approximate chi squared 0022 0023 % Copyright (C) Mike Brookes 1998 0024 % Version: $Id: v_randvec.m 10865 2018-09-21 17:22:45Z dmb $ 0025 % 0026 % VOICEBOX is a MATLAB toolbox for speech processing. 0027 % Home page: http://www.ee.ic.ac.uk/hp/staff/dmb/voicebox/voicebox.html 0028 % 0029 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0030 % This program is free software; you can redistribute it and/or modify 0031 % it under the terms of the GNU General Public License as published by 0032 % the Free Software Foundation; either version 2 of the License, or 0033 % (at your option) any later version. 0034 % 0035 % This program is distributed in the hope that it will be useful, 0036 % but WITHOUT ANY WARRANTY; without even the implied warranty of 0037 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 0038 % GNU General Public License for more details. 0039 % 0040 % You can obtain a copy of the GNU General Public License from 0041 % http://www.gnu.org/copyleft/gpl.html or by writing to 0042 % Free Software Foundation, Inc.,675 Mass Ave, Cambridge, MA 02139, USA. 0043 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0044 0045 % first sort out the input arguments 0046 0047 sm=size(m); 0048 if nargin<3 0049 c=ones(sm); % default to unit variance 0050 end 0051 sc=size(c); 0052 p=max(sm(2),sc(2)); % data dimension 0053 k=sm(1); % number of mixtures 0054 fullc=(length(sc)>2) || (sc(1)>k); 0055 if nargin<4 0056 mode='g'; % default to gaussian 0057 w=ones(k,1); 0058 else 0059 if ischar(w) 0060 mode = w; % w argument has been omitted 0061 w=ones(k,1); 0062 elseif nargin<5 0063 mode='g'; 0064 end 0065 end 0066 ty=mode(1); % ignore all but first character for type 0067 x=zeros(n,p); % initialize output array 0068 if sm(2)~=p 0069 m=repmat(m,1,p); % if m is given as a scalar 0070 end 0071 if sc(2) ~=p 0072 c=repmat(c,1,p); % if c is given as a scalar 0073 end 0074 q=sqrt(0.5); 0075 if k>1 0076 kx=v_randiscr(w,n); 0077 else 0078 kx=ones(n,1); 0079 end 0080 for kk=1:k 0081 nx=find(kx==kk); 0082 nn=length(nx); 0083 if nn % check if we need to generate any from mixture kk 0084 0085 % extract the mean and cov for this mixture 0086 0087 mm=m(kk,:); % mean vector 0088 if fullc % full covariance matrix 0089 cc=c(:,:,kk); 0090 if ty=='l' % lognormal distribution - convert mean and covariance 0091 cc=log(1+cc./(mm.'*mm)); 0092 mm=log(mm)-0.5*diag(cc).'; 0093 end 0094 else 0095 cc=c(kk,:); 0096 if ty=='l' % lognormal distribution - convert mean and covariance 0097 cc=log(1+cc(:).'./mm.^2); 0098 mm=log(mm)-0.5*cc; 0099 end 0100 end 0101 0102 % now generate nn complex or real values 0103 0104 if ty=='c' % generate complex values 0105 xx=q*randn(nn,p)+1i*q*randn(nn,p); % complex-valued unit variance values 0106 else 0107 xx=randn(nn,p); % real-valued unit variance values 0108 end; 0109 0110 % scale by the square root of the covariance matrix 0111 0112 if fullc % full covariance covariance 0113 [v,d]=eig((cc+cc')/2); % force covariance matrix to be hermitian 0114 xx=(xx.*repmat(sqrt(abs(diag(d))).',nn,1))*v'+repmat(mm,nn,1); % and also positive definite 0115 else 0116 xx=xx.*repmat(sqrt(abs(cc)),nn,1)+repmat(mm,nn,1); % different mean/cov for each column 0117 end 0118 x(nx,:)=xx; 0119 end 0120 end 0121 if ty=='l' % lognormal distribution 0122 x=exp(x); 0123 end 0124 if ~nargout 0125 if p==1 0126 if ty=='c' 0127 plot(real(x), imag(x),'+'); 0128 xlabel('Real'); 0129 ylabel('Imag'); 0130 else 0131 nbin=max(min(floor(sqrt(n)),50),5); 0132 hist(x,nbin); 0133 xlabel('X'); 0134 ylabel('Frequency'); 0135 end 0136 else 0137 [vv,iv]=sort(var(x,0,1)); 0138 iv=sort(iv(end-1:end)); 0139 plot(real(x(:,iv(1))), real(x(:,iv(2))),'+'); 0140 if ty=='c' 0141 xylab='Real[ x(%d) ]'; 0142 else 0143 xylab='x(%d)'; 0144 end 0145 xlabel(sprintf(xylab,iv(1))); 0146 ylabel(sprintf(xylab,iv(2))); 0147 end 0148 end 0149