Initial commit

This commit is contained in:
Sam Perry
2016-12-01 20:53:47 +00:00
commit 40bf0d1b60
6 changed files with 298 additions and 0 deletions
+43
View File
@@ -0,0 +1,43 @@
function [x_class, Perr] = classify()
phoneme_1 = load('MoGVarsK3.mat');
phoneme_2 = load('MoGVarsK3_p2.mat');
J = load('PB_data.mat');
% x = load('PB12.mat.mat');
x = [J.f1(J.phno == 1), J.f2(J.phno == 1)];
x1_min = min(J.f1)
x1_max = max(J.f1)
x2_min = min(J.f2)
x2_max = max(J.f2)
x = [linspace(x1_min, x1_max, 20)', linspace(x2_min, x2_max, 20)'];
[p1, d1] = liklihood(x, phoneme_1);
[p2, d2]= liklihood(x, phoneme_2);
x_class = mean(p1, 2) > mean(p2, 2)
y_class = mean(p1, 2) <= mean(p2, 2)
Perr = (0.5*sum(p1(x_class).*d1(x_class)))+(0.5*sum(p2(~x_class).*d2(~x_class)));
x_class + y_class'
function [Px, d] = liklihood(x, phoneme)
[n D] = size(x);
k = size(phoneme.mu, 2);
Px = zeros(n, k);
d = zeros(n, k);
p = phoneme.p;
mu = phoneme.mu;
s2 = phoneme.s2;
for i=1:k
Px(:,i) = p(i)*det(s2(:,:,i))^(-0.5)*exp(-0.5*sum((x'-repmat(mu(:,i),1,n))'*inv(s2(:,:,i)).*(x'-repmat(mu(:,i),1,n))',2));
d(:,i) = abs(sum((x'-repmat(mu(:,i),1,n))));
end
%Z = sum(Z, 2) / k;
%Px = mean(Px);
+105
View File
@@ -0,0 +1,105 @@
function demo_classification
% Generative classification demo: The Gaussian classifier with shared diagonal covariance.
% That is, a Gaussian Naive Bayes with shared covariance.
% Note, whenever the covariance is shared (i.e. the same for both classes) then the decision function is linear.
% This code is provided to complete Assignment 1 for the module Machine Learning (Extended)
% It (akong with Assignment 1) also serves as a starting point for you if you want to experiment with variations
% of the Gaussian Classifier such as non-diagonal class-covariance, non-shared class covariances etc.
% Generate some data in 2 classes
nInputs = 2; % dimensionality of data (number of attributes)
nClasses = 2; % number of classes
nSamples(1) = 200; nSamples(2) = 200; % number of points in each class
mu1 = [1;-1]; sigma1 = 1.2;
mu2 = [-1;1]; sigma2 = 1.8;
X = zeros(nInputs,max(nSamples),nClasses);
X(:,1:nSamples(1),1) = sigma1*randn(nInputs,nSamples(1)) + ...
repmat(mu1,1,nSamples(1)); % class 1 training data
X(:,1:nSamples(2),2) = sigma2*randn(nInputs,nSamples(2)) + ...
repmat(mu2,1,nSamples(2)); % class 2 training data
% Alternatively, load data
% targets
Y = [1,0]; % targets for class 1,2.
% Plot the data
clf, figure(1)
plot(X(1,1:nSamples(1),1),X(2,1:nSamples(1),1),'b+'); hold on % class 1
plot(X(1,1:nSamples(1),2),X(2,1:nSamples(1),2),'ro'); % class 2
title('2-CLASS DATA')
disp('Press a key to continue training the GAUSSIAN CLASSIFIER model ...'); pause
N = sum(nSamples); % total number of points
Sigma = zeros(nInputs);
for c = 1:nClasses
Nc = nSamples(c); Xc = X(:,1:Nc,c);
alpha(c) = Nc/N; % class prior
% estimate the mean of class c
m(:,c) = sum(Xc,2)/Nc; % Equivalent to: mean(Xc')
% estimate covariance of diagonal form for class c
% Sigma(:,:,c) = diag(diag((Xc-repmat(m(:,c),1,Nc))*(Xc-repmat(m(:,c),1,Nc))'))/Nc; % each class with its own covariance
% Equiv. to: cov(Xc')
Sigma = Sigma + (Xc-repmat(m(:,c),1,Nc))*(Xc-repmat(m(:,c),1,Nc))'; % instead, here we take a shared covariance
% this is estimated from all the points, not just from class c.
end%for c
% ML if all classes have same variance
Sigma = diag(diag(Sigma))/N; % diagonal approximation of the shared covariance estimate.
inv_Sigma = inv(Sigma);
figure(1)
plotgaus(m(:,1)',Sigma,'b');
plotgaus(m(:,2)',Sigma,'r');
% Computing & plotting the estimated decision boundary -- this is linear because the class covariances
% were taken to be equal in the two classes. In case you decide to experiment with each class having its own covariance then
% remove the following three lines.
delta_m = m(:,1)-m(:,2); mean_m = (m(:,1)+m(:,2))/2;
m
m(:,1)
m(:,2)
w = inv_Sigma*delta_m;
b = -delta_m'*inv_Sigma*mean_m + log(alpha(1)/(1-alpha(1)));
xmin = -3; xmax = 3;
hl2 = line([xmin,xmax],[(-b-w(1)*xmin)/w(2),(-b-w(1)*xmax)/w(2)]);
set(hl2,'color','g','linewidth',2);
drawnow
function [h] = plotgaus( mu, sigma, colspec );
% PLOTGAUS Plotting of a Gaussian contour
%
% PLOTGAUS(MU,SIGMA,COLSPEC) plots the mean and standard
% deviation ellipsis of the Gaussian process that has mean MU
% variance SIGMA, with color COLSPEC = [R,G,B].
%
% If you use only PLOTGAUS(MU,SIGMA), the default color is
% [0 1 1] (cyan).
%
if nargin < 3; colspec = [0 1 1]; end;
npts = 100;
stdev = sqrtm(sigma);
t = linspace(-pi, pi, npts);
t=t(:);
X = [cos(t) sin(t)] * stdev + repmat(mu,npts,1);
h(1) = line(X(:,1),X(:,2),'color',colspec,'linew',2);
h(2) = line(mu(1),mu(2),'marker','+','markersize',10,'color',colspec,'linew',2);
+17
View File
@@ -0,0 +1,17 @@
function main()
J = load('PB_data.mat');
% x = load('PB12.mat.mat');
x = [J.f1(J.phno == 2), J.f2(J.phno == 2)]
for i=1:10
[mu, s2, p] = mog(x, 6);
mu
s2
p
disp('Press a key to continue')
pause;
save('MoGVarsK6_p2.mat', 'mu', 's2', 'p')
end
% plot(J(:,1),J(:,2),'.')
View File
+60
View File
@@ -0,0 +1,60 @@
% Simple script to do EM for a mixture of Gaussians.
% -------------------------------------------------
% based on code from Rasmussen and Ghahramani
% (http://www.gatsby.ucl.ac.uk/~zoubin/course02/)
% Initialise parameters
function [mu, s2, p] = mog(x, k)
[n D] = size(x); % number of observations (n) and dimension (D)
p = ones(1,k)/k; % mixing proportions
mu = x(ceil(n.*rand(1,k)),:)'; % means picked randomly from data
s2 = zeros(D,D,k); % covariance matrices
niter=100; % number of iterations
% initialize covariances
for i=1:k
s2(:,:,i) = cov(x)./k; % initially set to fraction of data covariance
end
set(gcf,'Renderer','zbuffer');
clear Z;
try
% run EM for niter iterations
for t=1:niter,
fprintf('t=%d\r',t);
% Do the E-step:
for i=1:k
Z(:,i) = p(i)*det(s2(:,:,i))^(-0.5)*exp(-0.5*sum((x'-repmat(mu(:,i),1,n))'*inv(s2(:,:,i)).*(x'-repmat(mu(:,i),1,n))',2));
end
Z = Z./repmat(sum(Z,2),1,k);
% Do the M-step:
for i=1:k
mu(:,i) = (x'*Z(:,i))./sum(Z(:,i));
% We will fit Gaussians with diagonal covariances:
s2(:,:,i) = diag((x'-repmat(mu(:,i),1,n)).^2*Z(:,i)./sum(Z(:,i)));
% To fit general Gaussians use the line:
% s2(:,:,i) =
% (x'-repmat(mu(:,i),1,n))*(repmat(Z(:,i),1,D).*(x'-repmat(mu(:,i),1,n))')./sum(Z(:,i));
p(i) = mean(Z(:,i));
end
clf
hold on
plot(x(:,1),x(:,2),'.');
for i=1:k
plot_gaussian(2*s2(:,:,i),mu(:,i),i,11);
end
drawnow;
end
catch
disp('Numerical Error in Loop - Possibly Singular Matrix');
end;
+73
View File
@@ -0,0 +1,73 @@
% Plots a 2D or 3D 1 s.d. frame.
% 2D: 'n-1' divisions polar-wise,
% 3D: 'n-1' divisions each azimuthally and polar-wise,
% for a Gaussian with covariance 'covar' and mean 'mu',
% 'colour' can be any integer. M.Beal GMLC 26/03/99
%
% function plot_gaussian(covar,mu,col,n)
function [hh] = plot_gaussian(covar,mu,col,n)
if ~isempty(find(covar-covar == 0))
if size(mu,1) < size(mu,2), mu = mu'; end
if size(covar,1) == 3
theta = (0:1:n-1)'/(n-1)*pi;
phi = (0:1:n-1)/(n-1)*2*pi;
sx = sin(theta)*cos(phi);
sy = sin(theta)*sin(phi);
sz = cos(theta)*ones(1,n);
svect = [reshape(sx,1,n*n); reshape(sy,1,n*n); reshape(sz,1,n*n)];
epoints = sqrtm(covar) * svect + mu*ones(1,n*n);
ex = reshape(epoints(1,:),n,n);
ey = reshape(epoints(2,:),n,n);
ez = reshape(epoints(3,:),n,n);
colourset = [1 0 0; 0 1 0; 0 0 1; 1 1 0; 1 0 1; 0 1 1];
colour = colourset(mod(col-1,size(colourset,1))+1,:);
hh = mesh(ex,ey,ez, reshape(ones(n*n,1)*colour,n,n,3) );
hidden off
light
% for conversion to .eps figures, the 'mesh' command screws up the
% resolution. Therefore use these 4 lines insteasd of the previous 5.
% colourset = ['r'; 'g'; 'b'; 'y'; 'm'; 'c'];
% colour = colourset(mod(col-1,size(colourset,1))+1,:);
% plot3(epoints(1,:),epoints(2,:),epoints(3,:),colour)
% plot3(reshape(ex',1,n*n),reshape(ey',1,n*n),reshape(ez',1,n*n),colour)
else
theta = (0:1:n-1)/(n-1)*2*pi;
epoints = sqrtm(covar) * [cos(theta); sin(theta)] + mu*ones(1,n);
colourset = ['r'; 'g'; 'b'; 'y'; 'm'; 'c'];
colour = colourset(mod(col-1,size(colourset,1))+1,:);
hh = plot(epoints(1,:),epoints(2,:),colour,'LineWidth',3);
plot(mu(1,:),mu(2,:),[colour '.'])
end
else
fprintf('\nVery ill covariance matrix - not plotting this one\n')
end