Daniel@0: %DEMGAUSS Demonstrate sampling from Gaussian distributions.
Daniel@0: %
Daniel@0: %	Description
Daniel@0: %
Daniel@0: %	DEMGAUSS provides a simple illustration of the generation of data
Daniel@0: %	from Gaussian distributions. It first samples from a one-dimensional
Daniel@0: %	distribution using RANDN, and then plots a normalized histogram
Daniel@0: %	estimate of the distribution using HISTP together with the true
Daniel@0: %	density calculated using GAUSS.
Daniel@0: %
Daniel@0: %	DEMGAUSS then demonstrates sampling from a Gaussian distribution in
Daniel@0: %	two dimensions. It creates a mean vector and a covariance matrix, and
Daniel@0: %	then plots contours of constant density using the function GAUSS. A
Daniel@0: %	sample of points drawn from this distribution, obtained using the
Daniel@0: %	function GSAMP, is then superimposed on the contours.
Daniel@0: %
Daniel@0: %	See also
Daniel@0: %	GAUSS, GSAMP, HISTP
Daniel@0: %
Daniel@0: 
Daniel@0: %	Copyright (c) Ian T Nabney (1996-2001)
Daniel@0: 
Daniel@0: clc
Daniel@0: mean = 2; var = 5; nsamp = 3000;
Daniel@0: xmin = -10; xmax = 10; nbins = 30;
Daniel@0: disp('Demonstration of sampling from a uni-variate Gaussian with mean')
Daniel@0: dstring = [num2str(mean), ' and variance ', num2str(var), '.  ', ...
Daniel@0:     num2str(nsamp), ' samples are taken.'];
Daniel@0: disp(dstring);
Daniel@0: x = mean + sqrt(var)*randn(nsamp, 1);
Daniel@0: fh1 = figure;
Daniel@0: histp(x, xmin, xmax, nbins);
Daniel@0: hold on;
Daniel@0: axis([xmin xmax 0 0.2]);
Daniel@0: plotvals = linspace(xmin, xmax, 200)';
Daniel@0: probs = gauss(mean, var, plotvals);
Daniel@0: plot(plotvals, probs, '-r');
Daniel@0: xlabel('X')
Daniel@0: ylabel('Density')
Daniel@0: 
Daniel@0: disp(' ')
Daniel@0: disp('Press any key to continue')
Daniel@0: pause; 
Daniel@0: mu = [3 2];
Daniel@0: lam1 = 0.5;
Daniel@0: lam2 = 5.0;
Daniel@0: Sigma = lam1*[1,1]'*[1,1] + lam2*[1,-1]'*[1,-1];
Daniel@0: disp(' ')
Daniel@0: disp('Demonstration of sampling from a bi-variate Gaussian.  The mean is')
Daniel@0: dstring = ['[', num2str(mu(1)), ', ', num2str(mu(2)), ...
Daniel@0:       '] and the covariance matrix is'];
Daniel@0: disp(dstring)
Daniel@0: disp(Sigma);
Daniel@0: ngrid = 40;
Daniel@0: cmin = -5; cmax = 10; 
Daniel@0: cvals = linspace(cmin, cmax, ngrid);
Daniel@0: [X1, X2] = meshgrid(cvals, cvals);
Daniel@0: XX = [X1(:), X2(:)];
Daniel@0: probs = gauss(mu, Sigma, XX);
Daniel@0: probs = reshape(probs, ngrid, ngrid);
Daniel@0: 
Daniel@0: fh2 = figure;
Daniel@0: contour(X1, X2, probs, 'b');
Daniel@0: hold on
Daniel@0: 
Daniel@0: nsamp = 300;
Daniel@0: dstring = [num2str(nsamp), ' samples are generated.'];
Daniel@0: disp('The plot shows the sampled data points with a contour plot of their density.')
Daniel@0: samples = gsamp(mu, Sigma, nsamp);
Daniel@0: plot(samples(:,1), samples(:,2), 'or');
Daniel@0: xlabel('X1')
Daniel@0: ylabel('X2')
Daniel@0: grid off;
Daniel@0: 
Daniel@0: disp(' ')
Daniel@0: disp('Press any key to end')
Daniel@0: pause; 
Daniel@0: close(fh1);
Daniel@0: close(fh2);
Daniel@0: clear all;