%Does PCA and Fisher's Linear Discriminant (FLD) clear all; hold off; close all; N = 100; Dog = zeros(N, 2); %Make fake two class data Dog(:, 1) = 2*rand(N, 1) - 1; Dog(:, 2) = (Dog(:, 1) + 0.25*randn(N, 1)) + 0.5; Cat = zeros(N, 2); Cat(:, 1) = 1 - 2*rand(N, 1) + 0.5*randn(N, 1); Cat(:, 2) = (Cat(:, 1) - 0.4*rand(N, 1)); figure(1); plot(Dog(:, 1), Dog(:, 2), 'ro', Cat(:, 1), Cat(:, 2), 'bo'); hold on; axis equal; MD = mean(Dog); MC = mean(Cat); %Now do PCA on entire set of data Dta(1:N, :, :) = Dog; Dta((N+1):2*N, :, :) = Cat; CovBetween = cov(Dta); [PCA, EV] = eig(CovBetween); line(2*[-PCA(1, 2), PCA(1, 2)], 2*[-PCA(2, 2), PCA(2, 2)], 'Color', 'k', 'LineWidth', 2, 'LineStyle', '-'); %Now do Fisher Linear Discriminant CovDog = cov(Dog); CovCat = cov(Cat); CovWithin = (CovDog + CovCat)/2; [FLD, FV] = eig(inv(CovWithin)*CovBetween); NmFLD = FLD(:, 2)/norm(FLD(:, 2)); line([-NmFLD(1), NmFLD(1)], [-NmFLD(2), NmFLD(2)], 'Color', 'r', 'LineWidth', 2); FDog = Dog*NmFLD; %Project onto decision axis FCat = Cat*NmFLD; MnD = mean(FDog); MnC = mean(FCat); if MnD < MnC; %Make sure category #1 has biggest projection FDog = -FDog; FCat = - FCat; MnD = -MnD; MnC = -MnC; end; Correct = zeros(1, 100); for KK = 1:100; Decide(KK) = MnD + (MnC - MnD)*KK/100; % Candidate thresholds between MnD & MnC Correct(KK) = sum(FDog > Decide(KK)) + sum(FCat < Decide(KK)); end; [MX, II] = max(Correct); [MX, FF] = max(flip(Correct)); %get other end of plateau of maxima Theta = (Decide(II) + Decide(length(Correct) - FF))/2 Vthresh = [1; -NmFLD(1)/NmFLD(2)]; %Vector for threshold separatrix Tshift = Theta/NmFLD(2); Vthresh = Vthresh/norm(Vthresh); line(2*[-Vthresh(1), Vthresh(1)], 2*[-Vthresh(2) + Tshift, Vthresh(2) + Tshift], 'Color', 'r', 'LineWidth', 2, 'LineStyle', '--'); % XX = -2:0.1:2; % YY = -(NmFLD(1)*XX + NmFLD(3)*XX.^2)/NmFLD(2) - Theta/NmFLD(2); % plot(XX, YY, 'k', 'LineWidth', 2); NumCorrect = sum((FDog > Theta) + (FCat < Theta)) Errors = sum((FDog < Theta) + (FCat > Theta)) % figure(2); plot(Dog(:, 1), Dog(:, 2), 'ko', Cat(:, 1), Cat(:, 2), 'ko');