%% Linear Discriminant Analysis(LDA) DEMO         by      Po-Shen Kuo  11/04/2010
clc;close all;clear all;
%% Load Fisher's iris flower data set
load fisheriris;
%% Plot data 
figure(1);
gscatter(meas(:,1), meas(:,2), species,'rgb','osd');
xlabel('Sepal Length');ylabel('Sepal Width');
%% Caculate mean value
LDA_data=meas(:,:);  

data_mean=mean(LDA_data);  %  the mean value of whole data set 
data_mean=data_mean';

class_1_mean=mean(LDA_data(1:50,:));  %  the mean value of class 1 - setosa
class_1_mean=class_1_mean';

class_2_mean=mean(LDA_data(51:100,:));  %  the mean value of class 2 - versicolor
class_2_mean=class_2_mean';

class_3_mean=mean(LDA_data(101:150,:));  %  the mean value of class 3 - virginical
class_3_mean=class_3_mean';

mean=[class_1_mean,class_2_mean,class_3_mean,data_mean];
%% Determine the between-class scatter matrix(S_b)
N = 50;
S_b=0;

for i=1:3
    S_b = S_b + N * (mean(:,i)-mean(:,4)) * (mean(:,i)-mean(:,4))';
end

S_b
%% Determine the within-class scatter matrix(S_w)
class_1_cov=cov(LDA_data(1:50,:));  %  Find the covariance matrix of class 1
class_2_cov=cov(LDA_data(51:100,:));  %  Find the covariance matrix of class 2
class_3_cov=cov(LDA_data(101:150,:));  %  Find the covariance matrix of class 3

S_w = class_1_cov + class_2_cov + class_3_cov;
S_w
%% Find the eigenvector and eigenvalue of (S_w_inv * S_b)
S_w_inv=inv(S_w);
[v,d]=eig(S_w_inv*S_b);
v
d
%% Project data points 
y=v'*LDA_data';
y=y';
y
%% Plot
figure(2);gscatter(y(:,1),y(:,2),species,'rgb','osd')
figure(3);gscatter(y(:,1),zeros(150,1),species,'rgb','osd')  %  Project onto X-axis
