%projekt SNN 2012
%klasyfikator o 3 neuronach wyjściowych dla 3 klas
clear;
clc;

dane = load('SNN_A.TXT');
test = load('SNN_B.TXT');
tabela = zeros(10,7);
%--czy_trenowac=0 -> normalne działanie
%--czy_trenowac=1 -> wczytanie wcześniej zapisanej sieci
czy_trenowac=1;

if czy_trenowac==0;
load('siec.mat', 'siec');
end;

dane=dane';
test=test';
%--poniższą linię odkomentować by testować na zbiorze uczącym--
%test=dane;

wart_klasy=0.8;

% utworzenie klas
dane(3:5,:)=[wart_klasy*(dane(3,:)==1)+0.1;wart_klasy*(dane(3,:)==2)+0.1; wart_klasy*(dane(3,:)==3)+0.1];
test(3:5,:)=[wart_klasy*(test(3,:)==1)+0.1;wart_klasy*(test(3,:)==2)+0.1; wart_klasy*(test(3,:)==3)+0.1];

wart_klasy=wart_klasy+0.1;

margines_gorny=0.6*(wart_klasy);    
margines_dolny=0.4*(wart_klasy);  

indx_cl1=find(dane(3,:)>=margines_gorny);
indx_cl2=find(dane(4,:)>=margines_gorny);
indx_cl3=find(dane(5,:)>=margines_gorny);

%wizualizacja zbioru uczącego
% figure(1);
% plot(dane(1,indx_cl1),dane(2,indx_cl1),'ro');
% hold on;
% plot(dane(1,indx_cl2),dane(2,indx_cl2),'go');
% plot(dane(1,indx_cl3),dane(2,indx_cl3),'bo');
% hold off;
% legend("klasa 1","klasa 2","klasa 3");

%wartości min i max zbioru uczącego
we1_min=min(dane(1,:));
we1_max=max(dane(1,:));
we2_min=min(dane(2,:));
we2_max=max(dane(2,:));

for l=4:14 

if czy_trenowac==1;
%utworzenie sieci neuronowej
siec=newff([we1_min we1_max; we2_min we2_max], [l 3], {'logsig' 'logsig'}, 'traingd');
siec.trainParam.goal=0.01;
siec.trainParam.epochs=100;
%inicjalizacja wag wg Nguyena-Widrowa
initnw(siec, 1);
initnw(siec, 2);
%uczenie alg. BP
siec=train(siec, [dane(1,:); dane(2,:)], dane(3:5,:));
%uczenie alg. LM
siec.trainFcn = 'trainlm';
siec.trainParam.epochs=200;
siec=train(siec, [dane(1,:); dane(2,:)], dane(3:5,:));
save('siec.mat', 'siec');
end;

%indeksy dla zbioru testowego
dane_check=test(1:2,:);
indx_cl1_test=find(test(3,:)>=margines_gorny);
indx_cl2_test=find(test(4,:)>=margines_gorny);
indx_cl3_test=find(test(5,:)>=margines_gorny);

%wizualizacja zbioru testowego
% figure(2);
% plot(test(1,indx_cl1_test),test(2,indx_cl1_test),'ro');
% hold on;
% plot(test(1,indx_cl2_test),test(2,indx_cl2_test),'go');
% plot(test(1,indx_cl3_test),test(2,indx_cl3_test),'bo');
% hold off;
% legend("klasa 1","klasa 2","klasa 3");

wyjscie=sim(siec,dane_check);

%sprawdzanie przynależności do danej klasy
indx_cl1_net=find((wyjscie(1,:)>margines_gorny)&(wyjscie(2,:)<margines_dolny)&(wyjscie(3,:)<margines_dolny));
indx_cl2_net=find((wyjscie(1,:)<margines_dolny)&(wyjscie(2,:)>margines_gorny)&(wyjscie(3,:)<margines_dolny));
indx_cl3_net=find((wyjscie(1,:)<margines_dolny)&(wyjscie(2,:)<margines_dolny)&(wyjscie(3,:)>margines_gorny));

%odczyt wag sieci
wagi_n_ukrytych=siec.IW{1,1}
bias_n_ukrytych=siec.b{1}
wagi_n_wyjsciowych=siec.LW{2,1}
bias_n_wyjsciowych=siec.b{2}

%wartości min i max zbioru testowego
we1_min=min(test(1,:));
we1_max=max(test(1,:));
we2_min=min(test(2,:));
we2_max=max(test(2,:));

% %wizualizacja linii rozdzielających klasy
% figure(3);
% plot(test(1,indx_cl1_net),test(2,indx_cl1_net),'r.');
% hold on;
% plot(test(1,indx_cl2_net),test(2,indx_cl2_net),'g.');
% plot(test(1,indx_cl3_net),test(2,indx_cl3_net),'b.');
% hold off;
il_neuronow_ukrytych=siec.layers{1}.size;
% figure(3);
% hold on;
for n=1:il_neuronow_ukrytych;
  if abs(wagi_n_ukrytych(n,1))<abs(wagi_n_ukrytych(n,2))
    plot([we1_min we1_max],[-wagi_n_ukrytych(n,1)/wagi_n_ukrytych(n,2)*we1_min-bias_n_ukrytych(n)/wagi_n_ukrytych(n,2) -wagi_n_ukrytych(n,1)/wagi_n_ukrytych(n,2)*we1_max-bias_n_ukrytych(n)/wagi_n_ukrytych(n,2)],'k');
  else
    plot([-wagi_n_ukrytych(n,2)/wagi_n_ukrytych(n,1)*we2_min-bias_n_ukrytych(n)/wagi_n_ukrytych(n,1) -wagi_n_ukrytych(n,2)/wagi_n_ukrytych(n,1)*we2_max-bias_n_ukrytych(n)/wagi_n_ukrytych(n,1)],[we2_min we2_max],'k');
  end;
end;
hold off;
axis([we1_min we1_max we2_min we2_max]);

%transpozycja indeksów
indx_cl1=indx_cl1';
indx_cl2=indx_cl2';
indx_cl3=indx_cl3';
indx_cl_n1=indx_cl1_test';
indx_cl_n2=indx_cl2_test';
indx_cl_n3=indx_cl3_test';

%obliczenia true positive
TP=zeros(3,1);
TN=zeros(3,1);
T=zeros(3,1);  
FN=zeros(3,1);
FP=zeros(3,1);
F=zeros(3,1);  

il_punktow=length(test);

a = 1;
for i=1:length(indx_cl1_net)
    znaleziono = 0;
    for j=1:length(indx_cl_n1)
        if (indx_cl1_net(i) == indx_cl_n1(j))
            znaleziono = 1;
            break;
        end
    end
    if (znaleziono)
        TP(a)=TP(a)+1;
    end
end
for i=1:length(indx_cl_n1)
    znaleziono = 0;
    for j=1:length(indx_cl1_net)
        if (indx_cl_n1(i) == indx_cl1_net(j))
            znaleziono=1;
            break;
        end
    end
    if (znaleziono == 0)
        FN(a)=FN(a)+1;
    end
end


a = 2;
for i=1:length(indx_cl2_net)
    znaleziono = 0;
    for j=1:length(indx_cl_n2)
        if (indx_cl2_net(i) == indx_cl_n2(j))
            znaleziono = 1;
            break;
        end
    end
    if (znaleziono)
        TP(a)=TP(a)+1;
    end
end
for i=1:length(indx_cl_n2)
    znaleziono = 0;
    for j=1:length(indx_cl2_net)
        if (indx_cl_n2(i) == indx_cl2_net(j))
            znaleziono=1;
            break;
        end
    end
    if (znaleziono == 0)
        FN(a)=FN(a)+1;
    end
end
   
a = 3;
for i=1:length(indx_cl3_net)
    znaleziono = 0;
    for j=1:length(indx_cl_n3)
        if (indx_cl3_net(i) == indx_cl_n3(j))
            znaleziono = 1;
            break;
        end
    end
    if (znaleziono)
        TP(a)=TP(a)+1;
    end
end
for i=1:length(indx_cl_n3)
    znaleziono = 0;
    for j=1:length(indx_cl3_net)
        if (indx_cl_n3(i) == indx_cl3_net(j))
            znaleziono=1;
            break;
        end
    end
    if (znaleziono == 0)
        FN(a)=FN(a)+1;
    end
end

%punkty niezklasyfikowane
Nklas = zeros(3,1);
a = 1;
for i=1:length(indx_cl_n1)
    szuk = indx_cl_n1(i);
    znaleziono=[];
    znaleziono = find((indx_cl1_net(:)==szuk));
    if (isempty(znaleziono))
        znaleziono = find(indx_cl2_net(:)==szuk);
    end
    if (isempty(znaleziono))
        znaleziono = find(indx_cl3_net(:)==szuk);
    end
    if (isempty(znaleziono))
        Nklas(a)=Nklas(a)+1;
    end
end
a = 2;
for i=1:length(indx_cl_n2)
    szuk = indx_cl_n2(i);
    znaleziono=[];
    znaleziono = find((indx_cl1_net(:)==szuk));
    if (isempty(znaleziono))
        znaleziono = find(indx_cl2_net(:)==szuk);
    end
    if (isempty(znaleziono))
        znaleziono = find(indx_cl3_net(:)==szuk);
    end
    if (isempty(znaleziono))
        Nklas(a)=Nklas(a)+1;
    end
end
a = 3;
for i=1:length(indx_cl_n3)
    szuk = indx_cl_n3(i);
    znaleziono=[];
    znaleziono = find((indx_cl1_net(:)==szuk));
    if (isempty(znaleziono))
        znaleziono = find(indx_cl2_net(:)==szuk);
    end
    if (isempty(znaleziono))
        znaleziono = find(indx_cl3_net(:)==szuk);
    end
    if (isempty(znaleziono))
        Nklas(a)=Nklas(a)+1;
    end
end

pkt_klasy_1=length(indx_cl1_test);
pkt_klasy_2=length(indx_cl2_test);
pkt_klasy_3=length(indx_cl3_test);

%obliczanie false positive, true negative etc.
FP(1)=pkt_klasy_1-TP(1)-Nklas(1);
FP(2)=pkt_klasy_2-TP(2)-Nklas(2);
FP(3)=pkt_klasy_3-TP(3)-Nklas(3);

TN(1)=(pkt_klasy_2 + pkt_klasy_3) - FP(1);
TN(2)=(pkt_klasy_1 + pkt_klasy_3) - FP(2);
TN(3)=(pkt_klasy_1 + pkt_klasy_2) - FP(3);

FP_calosc = FP(1)+FP(2)+FP(3);
TP_calosc = TP(1)+TP(2)+TP(3);
FN_calosc = FN(1)+FN(2)+FN(3);
TN_calosc = TN(1)+TN(2)+TN(3);

T_proc=zeros(3,1);
F_proc=zeros(3,1);

Nklas_suma = 0;
for i=1:3
    Nklas_suma = Nklas_suma + Nklas(i);
end
Nklas_proc=zeros(3,1);

%czułość, specyficzność itp. - ocena szczegółowa
sensitivity=TP_calosc / (TP_calosc + FN_calosc) * 100
fp_rate=FP_calosc / (FP_calosc + TN_calosc) * 100
precision=TP_calosc / (TP_calosc + FP_calosc) * 100
accuracy=(TP_calosc + TN_calosc) / (TP_calosc + FN_calosc + FP_calosc + TN_calosc) * 100
specificity=TN_calosc / (FP_calosc + TN_calosc) * 100

T=TP;
F=FP;
for i=1:3
    T_proc(i) = T(i) / (T(i) + F(i) + Nklas(i))*100;
    F_proc(i) = F(i) / (T(i) + F(i) + Nklas(i))*100;
    Nklas_proc(i) = Nklas(i) / (T(i) + F(i) + Nklas(i))*100;
end

TP_calosc_proc=(TP_calosc) / il_punktow * 100;
FP_calosc_proc=(FP_calosc) / il_punktow * 100;
Nklas_suma_proc = Nklas_suma / il_punktow * 100;

j=l-3;
tabela(j,1)= j;
tabela(j,2)= TP_calosc_proc;
tabela(j,3)= TP_calosc;
tabela(j,4)= FP_calosc;
tabela(j,5)= FP_calosc_proc;
tabela(j,6)= Nklas_suma;
tabela(j,7)= Nklas_suma_proc;

TP_calosc
end

T
T_proc
F
F_proc
TP_calosc
% TP_calosc_proc=(TP_calosc) / il_punktow * 100
FP_calosc
% FP_calosc_proc=(FP_calosc) / il_punktow * 100
Nklas
Nklas_proc
Nklas_suma
% Nklas_suma_proc = Nklas_suma / il_punktow * 100

il_neuronow_ukrytych
il_punktow

tabela


% czerwona_1=size(indx_cl1_test)
% zielona_2=size(indx_cl2_test)
% niebieska_3=size(indx_cl3_test)