samplegaussian.m

內涵模糊理論與類神經網路的程式碼...提供初學者做研究參考

% ==========================================================
% 
%           Neural Networks A Classroom Approach
%                     Satish Kumar
%             Copyright Tata McGraw Hill, 2004
%
%   MATLAB code to generate probability contours and sample
%   a 2 - dimensional gaussian distribution
%                Reference: Table 7.2;Page 228
%
% ==========================================================

% ONLY one of these sets has been described in the text.

clear;
mu = [2 2]; %Set up the mean

figure;
subplot(3,2,1);

%Set up the covariance matrix and its inverse
alpha1= 0.5; 
alpha2= 4; 
covar = alpha1*[1,1]'*[1,1] + alpha2*[1,-1]'*[1,-1]
invcov = inv(covar);

gridres = 40; % Resolution of grid for evaluation
min = -5;       % Set up the limits for the mesh
max = 10;      % over which the contours will evaluate
data = linspace(min, max, gridres);

[X1, X2] = meshgrid(data, data);
X = [X1(:), X2(:)];  % Vector of x and y components from the grid
[n, d] = size(X);    % Get the number and dimension of the data

%Compute the probability distribution over the grid
%First compute the Mahalanobis distance
X = X - ones(n, 1)*mu;
exponent = sum(((X*invcov).*X), 2);
p = exp(-0.5*exponent)./sqrt((2*pi)^d*det(covar));
P = reshape(p, gridres, gridres);

contour(X1, X2, P, 'b');
hold on

%Now sample the function
numsamp = 300;

%Compute the eigenvectors and eigenvalues of the covariance matrix
[eta, lambda] = eig(covar);
coeffs = randn(numsamp, d)*sqrt(lambda); 
samples = coeffs*eta' + ones(numsamp, 1)*mu ; %Rotate and shif the samples

%Plot the samples
plot(samples(:,1), samples(:,2), 'k.');

xlabel('X1')
ylabel('X2')
axis square;
title('(a)');
subplot(3,2,2)
x1= reshape(X1, gridres, gridres);
x2= reshape(X2, gridres, gridres);
mesh(X1, X2, P);
axis([-5 10 -5 10 0 0.06]);
title('(b)');
%%%%%%%%%%%

subplot(3,2,3);

%Set up the covariance matrix and its inverse
alpha1= 4; 
alpha2= 0.5; 
covar = alpha1*[1,1]'*[1,1] + alpha2*[1,-1]'*[1,-1] 
invcov = inv(covar);

gridres = 40; % Resolution of grid for evaluation
min = -5;       % Set up the limits for the mesh
max = 10;      % over which the contours will evaluate
data = linspace(min, max, gridres);

[X1, X2] = meshgrid(data, data);
X = [X1(:), X2(:)];  % Vector of x and y components from the grid
[n, d] = size(X);    % Get the number and dimension of the data

%Compute the probability distribution over the grid
%First compute the Mahalanobis distance
X = X - ones(n, 1)*mu;
exponent = sum(((X*invcov).*X), 2);
p = exp(-0.5*exponent)./sqrt((2*pi)^d*det(covar));
P = reshape(p, gridres, gridres);

contour(X1, X2, P, 'b');
hold on

%Now sample the function
numsamp = 300;

%Compute the eigenvectors and eigenvalues of the covariance matrix
[eta, lambda] = eig(covar);
coeffs = randn(numsamp, d)*sqrt(lambda); 
samples = coeffs*eta' + ones(numsamp, 1)*mu ; %Rotate and shif the samples

%Plot the samples
plot(samples(:,1), samples(:,2), 'k.');

xlabel('X1')
ylabel('X2')
axis square;
title('(c)');

subplot(3,2,4)
x1= reshape(X1, gridres, gridres);
x2= reshape(X2, gridres, gridres);
mesh(X1, X2, P);
axis([-5 10 -5 10 0 0.06]);
title('(d)');
%%%%%%%%%%

subplot(3,2,5);

%Set up the covariance matrix and its inverse
alpha1= 4.0; 
alpha2= 4; 
covar = alpha1*[1,1]'*[1,1] + alpha2*[1,-1]'*[1,-1]
invcov = inv(covar);

gridres = 40; % Resolution of grid for evaluation
min = -5;       % Set up the limits for the mesh
max = 10;      % over which the contours will evaluate
data = linspace(min, max, gridres);

[X1, X2] = meshgrid(data, data);
X = [X1(:), X2(:)];  % Vector of x and y components from the grid
[n, d] = size(X);    % Get the number and dimension of the data

%Compute the probability distribution over the grid
%First compute the Mahalanobis distance
X = X - ones(n, 1)*mu;
exponent = sum(((X*invcov).*X), 2);
p = exp(-0.5*exponent)./sqrt((2*pi)^d*det(covar));
P = reshape(p, gridres, gridres);

contour(X1, X2, P, 'b');
hold on

%Now sample the function
numsamp = 300;

%Compute the eigenvectors and eigenvalues of the covariance matrix
[eta, lambda] = eig(covar);
coeffs = randn(numsamp, d)*sqrt(lambda); 
samples = coeffs*eta' + ones(numsamp, 1)*mu ; %Rotate and shif the samples

%Plot the samples
plot(samples(:,1), samples(:,2), 'k.');

xlabel('X1')
ylabel('X2')
axis square;
title('(e)');

subplot(3,2,6)
x1= reshape(X1, gridres, gridres);
x2= reshape(X2, gridres, gridres);
mesh(X1, X2, P);
axis([-5 10 -5 10 0 0.02]);
title('(f)');