% Klasifikace dat pomoci Bayesova klasifikatoru
% 
% Autori: RNDr. Eva Janousova, Mgr. & Mgr. Petr Dluhos 
% Masarykova univerzita, Brno
% 20/08/2014
%
% ----------------------------------------------------------------------- %


%% Datovy soubor a jeho vizualizace

% Datovy soubor:
X1=[[5,6,7]',[5,2,2]'] % matice pacientu
X2=[[9,7,8]',[4,3,2]'] % matice kontrolnich subjektu
x1_mean=mean(X1,1) % centroid pacientu
x2_mean=mean(X2,1) % centroid kontrol
x0=[7,4] % novy subjekt, ktery chceme klasifikovat

% Vizualizace datoveho souboru:
figure
plot(X1(:,1),X1(:,2),'r.','MarkerSize',20) % vykresleni pacientu
axis([4,10,1,6]) % nastaveni os
hold on
plot(X2(:,1),X2(:,2),'k.','MarkerSize',20) % vykresleni kontrol
plot(x1_mean(1),x1_mean(2),'r*','MarkerSize',10) % vykresleni centroidu pacientu
plot(x2_mean(1),x2_mean(2),'k*','MarkerSize',10) % vykresleni centroidu kontrol
plot(x0(1),x0(2),'b.','MarkerSize',20) % vykresleni noveho subjektu
hold off

% Vypocet kovariancnich matic:
S1=cov(X1) % kovariancni matice pacientu
S2=cov(X2) % kovariancni matice kontrol


%% Bayesuv klasifikator

% 1. vykresleni hustoty pro objem amygdaly:
figure, subplot(1,2,1)
plot(-5:0.1:10,0.5*normpdf(-5:0.1:10,x1_mean(2),sqrt(S1(2,2))),'r','LineWidth',2) % vykresleni hustoty pro skupinu pacientu, ktera je vynasobena apriorni psti 0.5
axis([-5,10,-0.01,0.20])
hold on
plot(-5:0.1:10,0.5*normpdf(-5:0.1:10,x2_mean(2),sqrt(S2(2,2))),'k','LineWidth',2) % vykresleni hustoty pro skupinu kontrol, ktera je vynasobena apriorni psti 0.5
plot(x0(2),0,'bo') % vykresleni noveho subjektu
plot([x0(2),x0(2)],[0,0.15],'b') % vykresleni, kam spada novy subjekt 
plot(x0(2),0.5*normpdf(x0(2),x1_mean(2),sqrt(S1(2,2))),'bo') % vypocitana hodnota podminene hustoty psti vyskytu bodu x0 ve tride pacientu
plot(x0(2),0.5*normpdf(x0(2),x2_mean(2),sqrt(S2(2,2))),'bo') % vypocitana hodnota podminene hustoty psti vyskytu bodu x0 ve tride kontrol
hold off
xlabel('Objem amygdaly')
% subjekt tedy bude zarazen do tridy kontrol

% co by se stalo, kdyby apriorni psti trid nebyly stejne:
subplot(1,2,2)
plot(-5:0.1:10,0.7*normpdf(-5:0.1:10,x1_mean(2),sqrt(S1(2,2))),'r','LineWidth',2) % vykresleni hustoty pro skupinu pacientu, ktera je vynasobena apriorni psti 0.5
axis([-5,10,-0.01,0.20])
hold on
plot(-5:0.1:10,0.3*normpdf(-5:0.1:10,x2_mean(2),sqrt(S2(2,2))),'k','LineWidth',2) % vykresleni hustoty pro skupinu kontrol, ktera je vynasobena apriorni psti 0.5
plot(x0(2),0,'bo') % vykresleni noveho subjektu
plot([x0(2),x0(2)],[0,0.15],'b') % vykresleni, kam spada novy subjekt 
plot(x0(2),0.7*normpdf(x0(2),x1_mean(2),sqrt(S1(2,2))),'bo') % vypocitana hodnota podminene hustoty psti vyskytu bodu x0 ve tride pacientu
plot(x0(2),0.3*normpdf(x0(2),x2_mean(2),sqrt(S2(2,2))),'bo') % vypocitana hodnota podminene hustoty psti vyskytu bodu x0 ve tride kontrol
hold off
xlabel('Objem amygdaly - jiné apriorní pravděpodobnosti')
% v tomto pripade bychom zaradili subjekt do tridy pacientu
% bud to takto vahujeme, nebo posouvame hranici

% 2. vykresleni hustoty pro objem hipokampu:
figure
plot(1:0.1:13,0.5*normpdf(1:0.1:13,x1_mean(1),sqrt(S1(1,1))),'r','LineWidth',2) % vykresleni hustoty pro skupinu pacientu, ktera je vynasobena apriorni psti 0.5
axis([1,13,-0.01,0.22])
hold on
plot(1:0.1:13,0.5*normpdf(1:0.1:13,x2_mean(1),sqrt(S2(1,1))),'k','LineWidth',2) % vykresleni hustoty pro skupinu kontrol, ktera je vynasobena apriorni psti 0.5
plot(x0(1),0,'bo') % vykresleni noveho subjektu
plot([x0(1),x0(1)],[0,0.22],'b') % vykresleni, kam spada novy subjekt 
plot(x0(1),0.5*normpdf(x0(1),x1_mean(1),sqrt(S1(1,1))),'bo') % vypocitana hodnota podminene hustoty psti vyskytu bodu x0 ve tride pacientu
plot(x0(1),0.5*normpdf(x0(1),x2_mean(1),sqrt(S2(1,1))),'bo') % vypocitana hodnota podminene hustoty psti vyskytu bodu x0 ve tride kontrol
hold off
xlabel('Objem hipokampu')
% nemuzeme rozhodnout, kam bude subjekt zarazen - nahodne ho zaradime do
% jedne z trid

% 3. vykresleni hustoty dvourozmerneho normalniho rozdeleni:
x1 = -3:.2:3; x2 = -3:.2:3; % zvoleni gridu hodnot
[x1d,x2d] = meshgrid(x1+x1_mean(1),x2+x1_mean(2)); % grid hodnot pro pacienty
Fd = mvnpdf([x1d(:) x2d(:)],x1_mean,S1); % hustota 2D norm. rozd. pro pacienty
Fd = reshape(Fd,length(x2),length(x1)); % preskalovani hustoty
[x1h,x2h] = meshgrid(x1+x2_mean(1),x2+x2_mean(2)); % grid hodnot pro kontroly
Fh = mvnpdf([x1h(:) x2h(:)],x2_mean,S2); % hustota 2D norm. rozd. pro kontroly
Fh = reshape(Fh,length(x2),length(x1)); % preskalovani hustoty
figure
surf(x1+x1_mean(1),x2+x1_mean(2),0.5*Fd,'FaceColor',[1 0 0]); % vykresleni hustoty pacientu
hold on
surf(x1+x2_mean(1),x2+x2_mean(2),0.5*Fh,'FaceColor',[0.5 0.5 0.5]); % vykresleni hustoty kontrol
plot3([x0(1) x0(1)],[x0(2) x0(2)],[0 0.05],'LineWidth',2)
xlabel('x1'); ylabel('x2'); 
hold off
% subjekt bude zarazen do tridy kontrol

