% IANTA14: Programa MATLAB para "Proyecto IANTA"
%
% Representación temporal y espectral del sonido de anuros
%
% 14: Clasificación con número variable de parámetros
%       El orden de los parámetros se calcula por separabilidad (IANTA13)
%

function IANTA14

clear all;
clear global all;
clc;

global PATRON;
global carpeta; %Carpeta de trabajo de los sonidos
global ParUtil; %Parámetros utilizados
    
%-------------- Variables generales
[FicheroPatron,roiPatron,FicheroResumenPatrones]=DefinePatron;
DefineVariablesGenerales;
carpeta='sonidos\Procesados\';
%carpeta='sonidos\kk\';
FicheroResumen=[carpeta,'ResumenParametrosVariables'];

OrdPar=[6,1,8,10,7,12,5,3,18,4,17,16,14,13,11,15,9,2];
[matpar,ClasePatronNum,NomPar]=GeneraDatosPatrones;

NumPar=5; %Número de parámetros utilizados
ParUtil=OrdPar(1:NumPar);
matpar=matpar(:,ParUtil);


%{
%-------------- Calcula la clase de sonido real de cada fichero
[NombreFichero,ClaseSonidoFichero]=ClasificaCarpetaporNombre;
cabecera={'Fichero','Real'};
nhoja=NumPar;
xlswrite(FicheroResumen,cabecera,nhoja);
xlswrite(FicheroResumen,NombreFichero',nhoja,'A2');
xlswrite(FicheroResumen,ClaseSonidoFichero',nhoja,'B2');
%}

%{
%-------------- Calcula la clase de sonido por ÁRBOL DE DECISIÓN
algoritmo=1;
cabecera={'Arbol Decision'};
display(cabecera);
clasificador=fitctree(matpar,ClasePatronNum);
[ClaseSonidoFichero,ProbClasFichero]=ClasificaCarpeta(algoritmo,clasificador);
EscribeColumnas(FicheroResumen,NumPar,algoritmo,...
    cabecera,ClaseSonidoFichero,ProbClasFichero);
%}

%{
%-------------- Calcula la clase de sonido por kNN (k vecinos más próximos)
algoritmo=2;
cabecera={'kNN'};
display(cabecera);
clasificador=fitcknn(matpar,ClasePatronNum);
[ClaseSonidoFichero,ProbClasFichero]=ClasificaCarpeta(algoritmo,clasificador);
EscribeColumnas(FicheroResumen,NumPar,algoritmo,...
    cabecera,ClaseSonidoFichero,ProbClasFichero);
%}

%{
%-------------- Calcula la clase de sonido por SVM (máquina de vectores soporte)
algoritmo=3;
cabecera={'SVM'};
display(cabecera);
clasificador=fitcMsvm(matpar,ClasePatronNum);
[ClaseSonidoFichero,ProbClasFichero]=ClasificaCarpeta(algoritmo,clasificador);
EscribeColumnas(FicheroResumen,NumPar,algoritmo,...
    cabecera,ClaseSonidoFichero,ProbClasFichero);
%}

%{
%-------------- Calcula la clase de sonido por REGRESIÓN LOGÍSTICA
algoritmo=4;
cabecera={'Regresión logística'};
display(cabecera);
clasificador=mnrfit(matpar,ClasePatronNum);
[ClaseSonidoFichero,ProbClasFichero]=ClasificaCarpeta(algoritmo,clasificador);
EscribeColumnas(FicheroResumen,NumPar,algoritmo,...
    cabecera,ClaseSonidoFichero,ProbClasFichero);
%}

%{
%-------------- Calcula la clase de sonido por redes neuronales
algoritmo=5;
cabecera={'Red neuronal'};
display(cabecera);
clasificador=feedforwardnet(10);
clasificador.trainParam.showWindow=0;
clasificador=train(clasificador,matpar',ClasePatronNum');
[ClaseSonidoFichero,ProbClasFichero]=ClasificaCarpeta(algoritmo,clasificador);
EscribeColumnas(FicheroResumen,NumPar,algoritmo,...
    cabecera,ClaseSonidoFichero,ProbClasFichero);
%}

%{
%-------------- Calcula la clase de sonido por función discriminante
algoritmo=6;
cabecera={'Función discriminante'};
display(cabecera);
clasificador=fitcdiscr(matpar,ClasePatronNum);
[ClaseSonidoFichero,ProbClasFichero]=ClasificaCarpeta(algoritmo,clasificador);
EscribeColumnas(FicheroResumen,NumPar,algoritmo,...
    cabecera,ClaseSonidoFichero,ProbClasFichero);
%}

%{
%-------------- Calcula la clase de sonido por redes bayesianas
algoritmo=7;
cabecera={'Red Bayesiana'};
display(cabecera);
clasificador=fitNaiveBayes(matpar,ClasePatronNum);
[ClaseSonidoFichero,ProbClasFichero]=ClasificaCarpeta(algoritmo,clasificador);
EscribeColumnas(FicheroResumen,NumPar,algoritmo,...
    cabecera,ClaseSonidoFichero,ProbClasFichero);
%}


%{
%-------------- Calcula la clase de sonido por distancia mínima
algoritmo=8;
cabecera={'Distancia mínima'};
display(cabecera);
clasificador=fitDistMin(matpar,ClasePatronNum);
[ClaseSonidoFichero,ProbClasFichero]=ClasificaCarpeta(algoritmo,clasificador);
EscribeColumnas(FicheroResumen,NumPar,algoritmo,...
    cabecera,ClaseSonidoFichero,ProbClasFichero);
%}

%{
%-------------- Calcula la clase de sonido por máxima verosimilitud
algoritmo=9;
cabecera={'Máxima Verosimilitud'};
display(cabecera);
clasificador=fitMaxVer(matpar,ClasePatronNum);
[ClaseSonidoFichero,ProbClasFichero]=ClasificaCarpeta(algoritmo,clasificador);
EscribeColumnas(FicheroResumen,NumPar,algoritmo,...
    cabecera,ClaseSonidoFichero,ProbClasFichero);
%}

%
%-------------- Gráficos algoritmos
nhoja=NumPar;
[~,~,hoja]=xlsread(FicheroResumen,nhoja);
tamhoja=size(hoja);
nfil=tamhoja(1);
ncol=tamhoja(2);
MatClas=cell2mat(hoja(2:nfil,2:ncol));
cabecera=hoja(1,:);
algoritmo=8;

%GraficaPerfilExito(MatClas(:,1),MatClas(:,2*algoritmo))
%GraficaProbExito(MatClas(:,1),MatClas(:,2*algoritmo),MatClas(:,2*algoritmo+1))
%GraficaCompAlg(MatClas,cabecera);
%GraficaCompNumPar(FicheroResumen);
%GraficaROCNumPar(FicheroResumen);
performance=EvaluacionClasificacion(hoja)

kk=1;
%}


%
% Función que genera los datos de los patrones
%
function [X,Y,NomPar]=GeneraDatosPatrones
% X: matriz de predictores (observación,parámetro)
% Y: vector de clases (observación)
% NomPar: nombre de paráemtros (parámetro)

global PATRON;

[FicheroPatron,roiPatron,FicheroResumenPatrones]=DefinePatron;
nesp=PATRON.nesp;
X=[];
ClasePatronTxt=[];
for iesp=1:nesp
    FicheroMuestras=FicheroPatron{iesp};
    nmuest=length(FicheroMuestras);
    for imuest=1:nmuest
        NombreFichero=FicheroMuestras{imuest};
        [~,~,hoja]=xlsread(NombreFichero);
        tamhoja=size(hoja);
        nfil=tamhoja(1);
        ncol=tamhoja(2);
        matparx=cell2mat(hoja(2:nfil,2:ncol-1));
        ClasePatronTxtX=hoja(2:nfil,ncol);
        X=[X;matparx];
        ClasePatronTxt=[ClasePatronTxt;ClasePatronTxtX];
    end %de una muestra
end % de una especie
Y=ClaseTxt2Num(ClasePatronTxt);
npar=(ncol-2)/2;
X=X(:,1:npar);
NomPar=hoja(1,2:npar+1);


%
% Función que clasifica los sonidos de un fichero 
%   (a partir del fichero Excel)
%
function [ClaseSonido,ProbClase]=ClasificaFichero...
    (NombreFichero,algoritmo,clasificador)
% NombreFichero: nombre del fichero a procesar (sin extensión)
% algoritmo: Identificador del algoritmo de clasificación
% clasificador: Identificador del clasificador utilizado
% ClaseSonido: Resultado de la clasificación
% ProbClase: Probabilidad de la clasificación (factor de mérito)

global PATRON;
global ParUtil; %Parámetros utilizados

nesp=PATRON.nesp;

[~,~,hoja]=xlsread(NombreFichero);
tamhoja=size(hoja);
nfil=tamhoja(1);
ncol=tamhoja(2);
matpar=cell2mat(hoja(2:nfil,2:ncol-1));
nobservaciones=size(matpar,1);

npar=ncol-2;
matpar=matpar(:,1:npar/2);
matpar=matpar(:,ParUtil);

switch algoritmo
    case 1 % Árbol de decisión
        [PredSonido,ProbSonido]=predict(clasificador,matpar);
        
    case 2 % kNN
        [PredSonido,ProbSonido]=predict(clasificador,matpar);

    case 3 % SVM
        [PredSonido,ProbSonido]=svmPredict(clasificador,matpar);

    case 4 % Regresión logística
        ProbSonido= mnrval(clasificador,matpar);
        [~,PredSonido]=max(ProbSonido,[],2);

    case 5 % Red neuronal
        PredSonido=round(clasificador(matpar')');
        PredSonido=min(PredSonido,nesp+1);
        PredSonido=max(PredSonido,1);
        ProbSonido=zeros(nobservaciones,nesp+1);
        for k=1:nobservaciones
            ProbSonido(k,PredSonido(k))=1;
        end
        
    case 6 % Función discriminante
        [PredSonido,ProbSonido]=predict(clasificador,matpar);

    case 7 % Red bayesiana
        [ProbSonido,PredSonido]=posterior(clasificador,matpar);

    case 8 % Distancia mínima
        [PredSonido,ProbSonido]=DistMinPredict(clasificador,matpar);

    case 9 % Máxima versoimilitud
        [PredSonido,ProbSonido]=MaxVerPredict(clasificador,matpar);
end

ContEsp=zeros(1,nesp);
for iesp=1:nesp
    ResComp=eq(PredSonido,iesp*ones(length(PredSonido),1));
    ContEsp(iesp)=sum(ResComp);
end
[~,ClaseSonido]=max(ContEsp);
ProbClase=ContEsp(ClaseSonido)/sum(ContEsp);
kk=1;


%
% Función que clasifica los sonidos de una carpeta 
%   (a partir de los nombres de los ficheros Excel)
%
function [NombreFichero,ClaseSonidoFichero]=ClasificaCarpetaporNombre
% NombreFichero: Nombre de cada fichero excel (vector)
% ClaseSonidoFichero: Clase de sonido de cada fichero excel (vector)

global carpeta; %Carpeta de trabajo de los sonidos

NombreFichero=ObtieneFicherosExcelCarpeta(carpeta);
nfich=length(NombreFichero);
for ifich=1:nfich
    nombre=NombreFichero{ifich};

    for i=length(nombre):-1:1
        car=nombre(i);
        if strcmp(car,'\')
            idSon=nombre(i+1:i+2);
            break;
        end
    end
    
    switch idSon
        case 'A1'
            ClaseSonidoFichero(ifich)=1;
        case 'A2'
            ClaseSonidoFichero(ifich)=1;
        case 'A3'
            ClaseSonidoFichero(ifich)=2;
        case 'B1'
            ClaseSonidoFichero(ifich)=3;
        case 'B2'
            ClaseSonidoFichero(ifich)=3;
    end

end %del fichero


%
% Función que clasifica los sonidos de una carpeta 
%   (a partir de los ficheros Excel)
%
function [ClaseSonidoFichero,ProbClasFichero]=...
    ClasificaCarpeta(algoritmo,clasificador)
% algoritmo: Identificador del algoritmo de clasificación
% clasificador: Identificador del clasificador utilizado
% ClaseSonidoFichero: Clase de sonido predicho de cada fichero excel (vector)
% ProbClasFichero: Probabilidad de la clase de sonido de cada fichero excel (vector)

global carpeta; %Carpeta de trabajo de los sonidos

NombreFichero=ObtieneFicherosExcelCarpeta(carpeta);
nfich=length(NombreFichero);
for ifich=1:nfich
    nombre=NombreFichero{ifich};
    disp(nombre);
    [ClaseSonidoFichero(ifich),ProbClasFichero(ifich)]=...
        ClasificaFichero(nombre,algoritmo,clasificador);
    disp(ClaseSonidoFichero(ifich));
end %del fichero


%
% Función que escribe las columnas del fichero excel correspondientes a un
%   algoritmo de clasificación
%
function EscribeColumnas(FicheroResumen,hoja,algoritmo,...
    cabecera,ClaseSonido,ProbClas)
% FicheroResumen: Nombre del fichero excel en el que se escribe
% hoja: Número de hoja en la que escribir
% algoritmo: Identificador del algoritmo de clasificación
% cabecera: Texo con la cabecera de las columnas
% ClaseSonido: Resultado de la clasificación
% ProbClas: Probabilidad de la clasificación (factor de mérito)

columna=2*algoritmo+1;
LetCol1=char('A'+columna-1);
LetCol2=char('A'+columna);
xlswrite(FicheroResumen,cabecera,hoja,[LetCol1,'1']);
xlswrite(FicheroResumen,ClaseSonido',hoja,[LetCol1,'2']);
xlswrite(FicheroResumen,ProbClas',hoja,[LetCol2,'2']);

%
% Función que traza el gráfico de perfíl de éxito de una técnica de
% clasificación
%
function GraficaPerfilExito(SonidoReal,ClaseSonido)
% ClaseSonido: Clase de sonido predicho de cada fichero excel (vector)
% SonidoReal: Clase de sonido real de cada fichero excel (vector)

global PATRON;

nesp=PATRON.nesp;
nexito=zeros(1,nesp);
nfichesp=zeros(1,nesp);
nfich=length(ClaseSonido);
for ifich=1:nfich
    nfichesp(SonidoReal(ifich))=nfichesp(SonidoReal(ifich))+1;
    if ClaseSonido(ifich)==SonidoReal(ifich)
        nexito(ClaseSonido(ifich))=nexito(ClaseSonido(ifich))+1;
    end
end

% Gráfico del perfil de éxito
theta=linspace(0,2*pi,nesp+1);
tasaexito=zeros(1,nesp+1);
tasaexito(1:nesp)=nexito./nfichesp;
tasaexito(nesp+1)=tasaexito(1);
rectangle('Position',[-1,-1,2,2],'Curvature',[1,1],'FaceColor','w');
axis off;
axis equal;
for tasa=[25,50,75]/100
    rectangle('Position',[-tasa,-tasa,2*tasa,2*tasa],'Curvature',[1,1],...
        'LineStyle',':');
end
for iesp=1:nesp
    th1=theta(iesp);
    xf=cos(th1);
    yf=sin(th1);
    color=PATRON.coloresp(iesp,:);
    line([0,xf],[0,yf],'Color',color);
    
    x1=tasaexito(iesp)*cos(th1);
    y1=tasaexito(iesp)*sin(th1);
    th2=theta(iesp+1);
    x2=tasaexito(iesp+1)*cos(th2);
    y2=tasaexito(iesp+1)*sin(th2);
    line([x1,x2],[y1,y2],'Color','k');
    
    delta=0.1;
    x1=x1-delta/2;
    y1=y1-delta/2;
    rectangle('Position',[x1,y1,delta,delta],'Curvature',[1,1],...
        'FaceColor',color);

    if xf>0
        ha='Left';
    elseif xf<0
        ha='Right';
    else
        ha='Center';
    end
    rt=1.1;
    xt=rt*cos(th1);
    yt=rt*sin(th1);
    label=[PATRON.NombreEspecie{iesp},' (',num2str(100*tasaexito(iesp),'%3.0f'),'%)'];
    text(xt,yt,label,'HorizontalAlignment',ha,'Color',color);
end
axis([-2,2,-1,1]);

exitoglobal=sum(nexito)/nfich*100;
title(['Exito global: ',num2str(exitoglobal,'%6.2f'),'%']);


%
% Función que traza el gráfico de la probabilidad de éxito 
%   de una técnica de clasificación
%
function GraficaProbExito(ClaseSonido,SonidoReal,ProbClase)
% ClaseSonido: Clase de sonido predicho de cada fichero excel (vector)
% SonidoReal: Clase de sonido real de cada fichero excel (vector)
% ProbClase: Probabilidad de la clase de sonido real de cada fichero excel (vector)

global PATRON;

nesp=PATRON.nesp;
nfichesp=zeros(1,nesp);
nfich=length(SonidoReal);
for ifich=1:nfich
    nfichesp(ClaseSonido(ifich))=nfichesp(ClaseSonido(ifich))+1;
end

% Dibujo de los áreas de cada especie real
xi=0;
for iesp=1:nesp
    if iesp~=nesp
        xf=sum(nfichesp(1:iesp))+0.5;
    else
        xf=sum(nfichesp(1:iesp))+1;
    end
    color=PATRON.coloresp(iesp,:);
    patch([xi,xf,xf,xi],[0,0,100,100],color,...
        'EdgeColor','none','FaceColor',color,'FaceAlpha',0.2);
    hold on;
    xi=xf;
end
axis([0,nfich+1,0,100]);

% Dibujo de la decisión de cada fichero
for ifich=1:nfich
    color=PATRON.coloresp(SonidoReal(ifich),:);
    plot([ifich,ifich],[0,100*ProbClase(ifich)],'Color',color);
    hold on;
    plot([ifich,ifich],100*[ProbClase(ifich),ProbClase(ifich)],'.','Color',color);
end %del fichero


%
% Función que traza el gráfico de comparación de número de parámetros
%
function GraficaCompNumPar(FicheroResumen)
% FicheroResumen: Nombre del fichero con el resumen de cálculos

global PATRON;
global ParUtil; %Parámetros utilizados

nesp=PATRON.nesp;

[~,~,hoja]=xlsread(FicheroResumen,1);
nfil=size(hoja,1);
ncol=size(hoja,2);
MatClas=cell2mat(hoja(2:nfil,2:ncol));
NumPar=length(ParUtil);
nfich=nfil-1;
nomAlg={'ArDec','kNN','SVM','RegLog','Neur','Discr','Bayes','DisMin','MaxVer'};
ordalg=[8,9,1:7];
nalg=length(ordalg);
ClaseReal=MatClas(:,1);
Merito=zeros(NumPar,nalg);
for inumpar=1:NumPar
    [~,~,hoja]=xlsread(FicheroResumen,inumpar);
    MatClas=cell2mat(hoja(2:nfil,2:ncol));
    nfichesp=zeros(nesp,1);
    nexito=zeros(nesp,nalg);
    %Cálculo del número de aciertos
    for ifich=1:nfich
        nfichesp(ClaseReal(ifich))=nfichesp(ClaseReal(ifich))+1;
        for ialg=1:nalg
            colalg=2*ialg;
            ClasePred=MatClas(:,colalg);
            if ClasePred(ifich)==ClaseReal(ifich)
                nexito(ClaseReal(ifich),ialg)=nexito(ClaseReal(ifich),ialg)+1;
            end
        end
    end
    %Cálculo de las tasas de éxito
    TasaExito=zeros(nesp,nalg);
    TasaExitoGlobal=zeros(1,nalg);
    RangoExitoGlobal=zeros(1,nalg);
    xText=cell(nalg,1);
    for ialg=1:nalg
        TasaExito(:,ialg)=100*nexito(:,ialg)./nfichesp;
        TasaExitoGlobal(ialg)=100*sum(nexito(:,ialg))./nfich;
        RangoExitoGlobal(ialg)=max(TasaExito(:,ialg))-min(TasaExito(:,ialg));
        x=100-TasaExitoGlobal(ialg);
        y=RangoExitoGlobal(ialg);
        Merito(inumpar,ialg)=sqrt(x^2+y^2);
    end
end %del número de parámetros

%
% Gráfico de la tasa de éxito de cada algoritmo
MarkerVec=['+','o','*','s','d','p'];
xText=cell(nalg,1);
for ix=1:nalg
    ialg=ordalg(ix);
    xText(ix)=nomAlg(ialg);
    imarker=rem(ix-1,length(MarkerVec))+1;
    plot(1:NumPar,Merito(:,ialg),...
        'Marker',MarkerVec(imarker));
    hold all;
end
legend(xText,'Location','BestOutside');

figure;
plot(1:NumPar,mean(Merito,2));
hold all;
plot(1:NumPar,min(Merito,[],2));
plot(1:NumPar,max(Merito,[],2));


%
% Función que traza el gráfico de comparación del ROC vs. el número de parámetros
%
function GraficaROCNumPar(FicheroResumen)
% FicheroResumen: Nombre del fichero con el resumen de cálculos

global PATRON;
global ParUtil; %Parámetros utilizados

NumPar=length(ParUtil);
nomAlg={'DisMin','MaxVer','ArDec','kNN','SVM','RegLog','Neur','Discr','Bayes'};
nalg=length(nomAlg);
ROC=zeros(NumPar,nalg);
F1=zeros(NumPar,nalg);
%Cálculo de la evaluación ROC
for inumpar=1:NumPar
    [~,~,hoja]=xlsread(FicheroResumen,inumpar);
    performance=EvaluacionClasificacion(hoja);
    PRC=performance(:,3);
    SNS=performance(:,4);
    SPC=1-performance(:,5);
    ROC(inumpar,:)=100*sqrt( SNS.^2 + SPC.^2 ) / sqrt(2);
    F1(inumpar,:)= 200* (PRC.*SNS) ./ (PRC+SNS); %F1 score
end %del número de parámetros

%
% Gráfico de la tasa de éxito de cada algoritmo
MarkerVec=['+','o','*','s','d','p'];
xText=cell(nalg,1);
%par = ROC;
par = F1;
for ialg=1:nalg
    xText(ialg)=nomAlg(ialg);
    imarker=rem(ialg-1,length(MarkerVec))+1;
    plot(1:NumPar,par(:,ialg),'Marker',MarkerVec(imarker));
    hold all;
end
xlabel('Number of parameters');
%ylabel('ROC (%)');
ylabel('F_1 (%)');
legend(xText,'Location','BestOutside');
ylim([0,100]);

figure;
%par = ROC;
par = F1;
plot(1:NumPar,max(par,[],2));
hold all;
plot(1:NumPar,mean(par,2));
algopt=3; %Algoritmo óptimo 
plot(1:NumPar,par(:,algopt));
plot(1:NumPar,min(par,[],2));
xlabel('Number of parameters');
%ylabel('ROC (%)');
ylabel('F_1 (%)');
leyenda=[{'Max'},{'Mean'},{'Best'},{'Min'}];
legend(leyenda,'Location','SouthEast');
ylim([0,100]);

%
% Función que calcula la evaluación de las técnicas de clasificación
%
function performance=EvaluacionClasificacion(hoja)
% hoja: Hoja con datos de clase de sonido para cada fichero excel y cada algoritmo
% performance: matriz de prestaciones (algoritmo,indicador)
%               1: Accuracy
%               2: Error Rate
%               3: Precision
%               4: Sensitivity
%               5: Specificity

tamhoja=size(hoja);
nfil=tamhoja(1);
ncol=tamhoja(2);
MatClas=cell2mat(hoja(2:nfil,2:ncol));
SonidoReal=MatClas(:,1);

%nalg=(ncol-2)/2;
ordalg=[8,9,1:7];
nalg=length(ordalg);
performance=zeros(nalg,7);
for ix=1:nalg
    ialg=ordalg(ix);
    ClaseSonido=MatClas(:,2*ialg);

    % Cálculo de la matriz de confusión (valores absolutos)
    confusion=confusionmat(SonidoReal,ClaseSonido);
    P=sum(confusion,2); %Positive
    TP=diag(confusion); %True Positive
    FN=P-TP; %False Negative

    PP=sum(confusion,1)'; %Predicted Positive
    FP=PP-TP; %False Positive
    N=sum(P)-P; %Negative
    TN=N-FP; %True Negative

    % Prestaciones para cada clase
    ACCx=(TP+TN)./(TP+TN+FP+FN); %Accuracy
    ERRx=(FP+FN)./(TP+TN+FP+FN); %Error Rate
    PRCx=TP./(TP+FP); %Precision
    SNSx=TP./(TP+FN); %Sensitivity
    SPCx=TN./(FP+TN); %Specificity
    PRCx(isnan(PRCx)) = 0; %Elimina NaN

    % Prestaciones para el conjunto de las clases
    ACCm = mean(ACCx);
    ERRm = mean(ERRx);
    PRCm = mean(PRCx);
    SNSm = mean(SNSx);
    SPCm = mean(SPCx);
    ROCm = sqrt( SNSm.^2 + SPCm.^2 ) / sqrt(2); %ROC
    F1m = 2* (PRCm.*SNSm) ./ (PRCm+SNSm); %F1 score

    performance(ix,1)=ACCm; %Accuracy
    performance(ix,2)=ERRm; %Error Rate
    performance(ix,3)=PRCm; %Precision
    performance(ix,4)=SNSm; %Sensitivity
    performance(ix,5)=SPCm; %Specificity
    performance(ix,6)=ROCm; %ROC
    performance(ix,7)=F1m; %F1 score
end
%display(performance);


%
% Función que traza el gráfico de comparación de algoritmos
%
function GraficaCompAlg(hoja,cabecera)
% hoja: Hoja con datos de clase de sonido para cada fichero excel y cada algoritmo
% cabecera: Cabecera de las columnas de las hojas (vector)

global PATRON;

nesp=PATRON.nesp;
tamhoja=size(hoja);
nfil=tamhoja(1);
ncol=tamhoja(2);
nfich=nfil;
nomAlg={'ArDec','kNN','SVM','RegLog','Neur','Discr','Bayes','DisMin','MaxVer'};
ordalg=[8,9,1:7];
nalg=length(ordalg);
ClaseReal=hoja(:,1);
nfichesp=zeros(nesp,1);
nexito=zeros(nesp,nalg);
for ifich=1:nfich
    nfichesp(ClaseReal(ifich))=nfichesp(ClaseReal(ifich))+1;
    for ix=1:nalg
        ialg=ordalg(ix);
        colalg=2*ialg;
        ClasePred=hoja(:,colalg);
        if ClasePred(ifich)==ClaseReal(ifich)
            nexito(ClaseReal(ifich),ix)=nexito(ClaseReal(ifich),ix)+1;
        end
    end
end

TasaExito=zeros(nesp,nalg);
TasaExitoGlobal=zeros(1,nalg);
RangoExitoGlobal=zeros(1,nalg);
xText=cell(nalg,1);
for ix=1:nalg
    TasaExito(:,ix)=100*nexito(:,ix)./nfichesp;
    TasaExitoGlobal(ix)=100*sum(nexito(:,ix))./nfich;
    RangoExitoGlobal(ix)=max(TasaExito(:,ix))-min(TasaExito(:,ix));
    xText{ix}=nomAlg{ordalg(ix)};
end

%
% Gráfico de la tasa de éxito de cada algoritmo
bar(1:nalg,TasaExitoGlobal,'c');
hold all;
for iesp=1:nesp
    plot(1:nalg,TasaExito(iesp,:),'*');
end
axis([0,nalg+1,0,110]);
set(gca,'XTickLabel',xText);
ylabel('Tasa de éxito');
%}

% Gráfico de la tasa de éxito vs. Rango de éxito
figure;
MarkerVec=['+','o','*','s','d','p'];
for ialg=1:nalg
    imarker=rem(ialg-1,length(MarkerVec))+1;
    plot(100-TasaExitoGlobal(ialg),RangoExitoGlobal(ialg),...
        'Marker',MarkerVec(imarker),'LineStyle','none');
    hold all;
end
legend(xText,'Location','Best');
xlabel('Tasa de error');
ylabel('Rango de tasa de error');
axis([0,100,0,100]);


%
% Obtiene los nombres de los ficheros Excel de una carpeta
%
function NombreFichero=ObtieneFicherosExcelCarpeta(carpeta)
% carpeta: nombre de la carpeta a analizar
% NombreFichero: nombre de los ficheros (vector de celdas)

dirant=cd;
cd(carpeta);
directorio=dir;
cd(dirant);

nficheros=length(directorio);
for i=1:nficheros
    nombres{i}=directorio(i).name;
end

[NombresOrd,iNombresOrd]=sort(nombres);
nfich=0;
for i=3:nficheros % Evita el '.' y el '..'
    j=iNombresOrd(i);
    NombreCompleto=directorio(j).name;
    longnom=length(NombreCompleto);
    idSon=NombreCompleto(1:3);
    switch idSon
        case 'A1-'
            idSonCorrecto=1;
        case 'A2-'
            idSonCorrecto=1;
        case 'A3-'
            idSonCorrecto=1;
        case 'B1-'
            idSonCorrecto=1;
        case 'B2-'
            idSonCorrecto=1;
        otherwise
            idSonCorrecto=0;
    end
    for k=longnom:-1:1
        if NombreCompleto(k)=='.'
            extension=NombreCompleto(k+1:longnom);
            break;
        end
    end
    if idSonCorrecto && strcmp(extension,'xls')
        nfich=nfich+1;
        NombreFichero{nfich}=[carpeta,NombreCompleto(1:k-1)];
    end
end


%
% Función que define las variables generales del programa
%
function DefineVariablesGenerales

global GENERAL; % Estructura de datos generales
    %GENERAL.hopSize: Duración de un "salto" de tiempo en el espectrograma (MP7)
    %GENERAL.msegSize: Duración de un microsegmento
    %GENERAL.frPrelev: Frecuencias de la banda relevante de potencias [baja,alta]
    %GENERAL.nForm: Número máximo de formantes considerados
    %GENERAL.NomParFrame: Nombre de los parámetros de un frame
    %GENERAL.NomCorParFrame: Nombre corte de los parámetros de un frame
    %GENERAL.ParSigEsp: Nombre de los parámetros significativos para determinar una especie en un frame
    %GENERAL.NumParPrimarios: Número de parámetros de frame primarios
    %GENERAL.NumParSecundarios: Número de parámetros de frame secundarios
    %GENERAL.NumParTerciarios: Número de parámetros de frame terciarios
    %GENERAL.gradoLPC: grado del polinomio usado en el algoritmo LPC
    %GENERAL.frFilt: Frecuencias de corte del filtr de preprocesamiento[baja,alta]
global PATRON;
    

GENERAL.hopSize=0.01; %Duración de un "salto" de tiempo en el espectrograma (MP7)
GENERAL.msegSize=0.1; %Duración de un microsegmento
GENERAL.frPrelev=[500,5000]; %Frecuencias de la banda relevante de potencias
GENERAL.nForm=3; %Número de formantes
GENERAL.gradoLPC=30; %Grado del polinomio usado en el algoritmo LPC
GENERAL.frFilt=[300,10000]; %Frecuencias de corte del filtro de preprocesamiento[baja,alta]

nesp=PATRON.nesp;
nform=GENERAL.nForm;
i=0;
i=i+1; GENERAL.NomParFrame{i}='Pt';
i=i+1; GENERAL.NomParFrame{i}='PRelev';
i=i+1; GENERAL.NomParFrame{i}='CP';
i=i+1; GENERAL.NomParFrame{i}='DE';
i=i+1; GENERAL.NomParFrame{i}='Pl';
i=i+1; GENERAL.NomParFrame{i}='Pitch';
i=i+1; GENERAL.NomParFrame{i}='Ra';
i=i+1; GENERAL.NomParFrame{i}='Fla';
for k=1:nform
    i=i+1; GENERAL.NomParFrame{i}=['FreqForm',num2str(k)];
    i=i+1; GENERAL.NomParFrame{i}=['AbForm',num2str(k)];
end
i=i+1; GENERAL.NomParFrame{i}='CA';
i=i+1; GENERAL.NomParFrame{i}='DeA';
i=i+1; GENERAL.NomParFrame{i}='DiA';
i=i+1; GENERAL.NomParFrame{i}='VA';
GENERAL.NumParPrimarios=i;

i=i+1; GENERAL.NomParFrame{i}='PtSigma';
i=i+1; GENERAL.NomParFrame{i}='PRelevSigma';
i=i+1; GENERAL.NomParFrame{i}='CPSigma';
i=i+1; GENERAL.NomParFrame{i}='DESigma';
i=i+1; GENERAL.NomParFrame{i}='PlSigma';
i=i+1; GENERAL.NomParFrame{i}='PitchSigma';
i=i+1; GENERAL.NomParFrame{i}='RaSigma';
i=i+1; GENERAL.NomParFrame{i}='FlaSigma';
for k=1:nform
    i=i+1; GENERAL.NomParFrame{i}=['FreqFormSigma',num2str(k)];
    i=i+1; GENERAL.NomParFrame{i}=['AbFormSigma',num2str(k)];
end
i=i+1; GENERAL.NomParFrame{i}='CASigma';
i=i+1; GENERAL.NomParFrame{i}='DeASigma';
i=i+1; GENERAL.NomParFrame{i}='DiASigma';
i=i+1; GENERAL.NomParFrame{i}='VASigma';
GENERAL.NumParSecundarios=i-GENERAL.NumParPrimarios;

for k=1:nesp
    i=i+1; GENERAL.NomParFrame{i}=['VrEsp',num2str(k)];
    i=i+1; GENERAL.NomParFrame{i}=['DisEsp',num2str(k)];
end
i=i+1; GENERAL.NomParFrame{i}='Relev';
GENERAL.NumParTerciarios=i-GENERAL.NumParPrimarios-GENERAL.NumParSecundarios;

i=0;
i=i+1; GENERAL.NomCorParFrame{i}='Pt';
i=i+1; GENERAL.NomCorParFrame{i}='PR';
i=i+1; GENERAL.NomCorParFrame{i}='CP';
i=i+1; GENERAL.NomCorParFrame{i}='DE';
i=i+1; GENERAL.NomCorParFrame{i}='Pl';
i=i+1; GENERAL.NomCorParFrame{i}='Pitch';
i=i+1; GENERAL.NomCorParFrame{i}='Ra';
i=i+1; GENERAL.NomCorParFrame{i}='Fla';
for k=1:nform
    i=i+1; GENERAL.NomCorParFrame{i}=['FrF',num2str(k)];
    i=i+1; GENERAL.NomCorParFrame{i}=['AbF',num2str(k)];
end
i=i+1; GENERAL.NomCorParFrame{i}='CA';
i=i+1; GENERAL.NomCorParFrame{i}='DeA';
i=i+1; GENERAL.NomCorParFrame{i}='DiA';
i=i+1; GENERAL.NomCorParFrame{i}='VA';

i=i+1; GENERAL.NomCorParFrame{i}='Pt~';
i=i+1; GENERAL.NomCorParFrame{i}='PR~';
i=i+1; GENERAL.NomCorParFrame{i}='CP~';
i=i+1; GENERAL.NomCorParFrame{i}='DE~';
i=i+1; GENERAL.NomCorParFrame{i}='Pl~';
i=i+1; GENERAL.NomCorParFrame{i}='Pitch~';
i=i+1; GENERAL.NomCorParFrame{i}='Ra~';
i=i+1; GENERAL.NomCorParFrame{i}='Fla~';
for k=1:nform
    i=i+1; GENERAL.NomCorParFrame{i}=['FF~',num2str(k)];
    i=i+1; GENERAL.NomCorParFrame{i}=['AB~',num2str(k)];
end
i=i+1; GENERAL.NomCorParFrame{i}='CA~';
i=i+1; GENERAL.NomCorParFrame{i}='DeA~';
i=i+1; GENERAL.NomCorParFrame{i}='DiA~';
i=i+1; GENERAL.NomCorParFrame{i}='VA~';

for k=1:nesp
    i=i+1; GENERAL.NomCorParFrame{i}=['VrE',num2str(k)];
    i=i+1; GENERAL.NomCorParFrame{i}=['DiE',num2str(k)];
end
i=i+1; GENERAL.NomCorParFrame{i}='Relev';


%
% Función que define los patrones
%
function [FicheroPatron,roiPatron,FicheroResumenPatrones]=DefinePatron
% FicheroPatron: Nombre de los ficheros Patrón (celda de mxn)
% roiPatron: Regiones de Interés (ROI) de los ficheros Patrón (celda de mxn)
% FicheroResumenPatrones: Nombre del fichero resumen de patrones (global PATRON)

global PATRON; % Estructura de datos de patrones
% Datos de cada patron de especie (se aproxima por una mezcla de N gaussianas)
    %PATRON.nesp: Número de especies en los patrones
    %PATRON.NombreEspecie: Nombre de las especies en los patrones
    %PATRON.coloresp: Color usado para las especies (matriz [tipo región,especie])
    %PATRON.ng: Número de gaussianas en la mezcla de funciones de distribución
    %PATRON.nt: Número de frames en los patrones (matriz [especie,muestra])
   

PATRON.nesp=3;
PATRON.NombreEspecie{1}='Sapo corredor';
PATRON.NombreEspecie{2}='Sapo corredor: suelta';
PATRON.NombreEspecie{3}='Sapo partero común';
PATRON.NombreEspecieLatin{1}='Epidalea Calamita';
PATRON.NombreEspecieLatin{2}='Epidalea Calamita: suelta';
PATRON.NombreEspecieLatin{3}='Alytes obstetricans';
PATRON.ng=2;
tablaColores=[0.00,0.00,1.00;...
              0.00,0.50,0.00;...
              1.00,0.00,0.00;...
              0.75,0.00,0.75];
PATRON.coloresp=tablaColores(1:PATRON.nesp+1,:);
FicheroResumenPatrones='sonidos\Procesados\Patrones\IANTA-PATRON';
FicheroPatron=cell(1,PATRON.nesp);
FicheroPatron{1}={'sonidos\Procesados\Patrones\0501_0-10',...
                  'sonidos\Procesados\Patrones\0707_0-10'};
FicheroPatron{2}={'sonidos\Procesados\Patrones\0503'};
FicheroPatron{3}={'sonidos\Procesados\Patrones\1325',...
                  'sonidos\Procesados\Patrones\1330'};

roiPatron=cell(1,PATRON.nesp);
roiPatron{1}={[  21,  78;
                127, 231;
                264, 366;
                395, 489;
                523, 628;
                664, 759;
                798, 900;
                933,1035],...
              [  20,  44;
                 73, 127;
                150, 194;
                221, 270;
                294, 339;
                362, 411;
                435, 484;
                507, 545;
                569, 611;
                635, 679;
                704, 747;
                771, 821;
                845, 887;
                910, 950;
                976,1018]};

roiPatron{2}={[  42,  58;
                128, 136;
                377, 388;
                745, 753;
                795, 812;
               1019,1034;
               1345,1356;
               1403,1418;
               1463,1479;
               1570,1577;
               1632,1642;
               1858,1865;
               1915,1925;
               2241,2251;
               2303,2313;
               2366,2380;
               2518,2531;
               2574,2582;
               2810,2817;
               2858,2873]};
              
roiPatron{3}={[  73,  89;
                150, 164;
                497, 513;
                876, 892;
               1433,1450;
               1808,1825;
               2206,2222;
               2556,2572;
               2933,2948;
               3406,3422;
               3677,3699;
               3721,3737],...
              [ 182, 194;
                263, 277;
                696, 711;
               1226,1239;
               1737,1751;
               2374,2388;
               2919,2934;
               3436,3451;
               4001,4015;
               4487,4501;
               4990,5005]};


%
% Función que convierte de la clase de un sonido de texto a número
%
function ClaseNum=ClaseTxt2Num(ClaseTxt)
% ClaseTxt: Clase de sonido en formato texto (vector)
% ClaseNum: Clase de sonido en formato número (vector)

global PATRON;
nesp=PATRON.nesp;
nelem=length(ClaseTxt);
ClaseNum=(nesp+1)*ones(nelem,1);
for k=1:nelem
    for iesp=1:nesp
        if strcmp(ClaseTxt(k),PATRON.NombreEspecie{iesp});
            ClaseNum(k)=iesp;
        end
    end
end


%
% Función que obtinen modelo SVM con múltiples clases
%
function clasificador=fitcMsvm(X,Y)
% X: matriz de predictores (observación,parámetro)
% Y: vector de clases (observación)
% clasificador: clasificador; vector (especies)

global PATRON;
nesp=PATRON.nesp;
nobservaciones=length(Y);
Y2=zeros(nobservaciones,1);
clasificador=cell(nesp+1,1);
for iesp=1:nesp+1
    for k=1:nobservaciones
        if Y(k)==iesp
            Y2(k)='A';
        else
            Y2(k)='B';
        end
    end %de las observaciones
    SVMModel=fitcsvm(X,Y2,'Standardize',true);
    CompactSVMModel=compact(SVMModel);
    clasificador{iesp}=fitPosterior(CompactSVMModel,X,Y2);
end %de las especies
save('SVM','clasificador');
kk=1;


%
% Función que clasifica mediante SVM con múltiples clases
%
function [Z,ProbZ]=svmPredict(clasificador,X)
% X: matriz de predictores (observación,parámetro)
% clasificador: clasificador; vector (clases)
% Z: vector de predicciones (observación)
% ProbZ: Probabilidad de las predicciones: matriz (observación,especie)

global PATRON;
nesp=PATRON.nesp;
nobservaciones=length(X);
ProbZ2=zeros(nobservaciones,2,nesp+1);
for iesp=1:nesp+1
    modeloSVM=clasificador{iesp};
    [~,ProbZ2(:,:,iesp)]=predict(modeloSVM,X);
end

ProbZ(:,:)=ProbZ2(:,1,:);
Z=zeros(nobservaciones,1);
for k=1:nobservaciones
    ProbZ(k,:)=ProbZ(k,:)/sum(ProbZ(k,:));
    [~,Z(k)]=max(ProbZ(k,:));
end
kk=1;


%
% Función que obtiene el clasificador por máxima verosimilitud
%
function clasificador=fitMaxVer(X,Y)
% X: matriz de predictores (observación,parámetro)
% Y: vector de clases (observación)
% clasificador: clasificador; vector (clases)

global PATRON;
nesp=PATRON.nesp;
ng=PATRON.ng;
npar=size(X,2);
gmd=cell(nesp+1);
for iesp=1:nesp+1
    Xesp=X(Y==iesp,:);
    GMDini.PComponents=0.5*ones(1,ng);
    medini=mean(Xesp);
    varini=var(Xesp);
    for i=1:npar
        GMDini.mu(:,i)=linspace(medini(i)-sqrt(varini(i)),medini(i)+sqrt(varini(i)),ng)';
    end
    for i=1:ng
        GMDini.Sigma(:,:,i)=cov(Xesp);
    end
    warning off all;
    %gmd{iesp}=fitgmdist(param,ng,'Start',GMDini);
    try
        gmd{iesp}=fitgmdist(Xesp,ng,'Start',GMDini);
    catch exception
        gmd{iesp}=gmdistribution(GMDini.mu,GMDini.Sigma,GMDini.PComponents);
    end
    warning on all;

end
clasificador=gmd;


%
% Función que clasifica mediante Máxima Verosimilitud
%
function [Z,ProbZ]=MaxVerPredict(clasificador,X)
% X: matriz de predictores (observación,parámetro)
% clasificador: clasificador
% Z: vector de predicciones (observación)
% ProbZ: Probabilidad de las predicciones: matriz (observación,especie)

global PATRON;
nesp=PATRON.nesp;
ng=PATRON.ng;
nobservaciones=size(X,1);
gmdEsp=clasificador;
Z=zeros(nobservaciones,1);
ProbZ=zeros(nobservaciones,nesp+1);
verosimilitud=zeros(nesp+1,1);

for k=1:nobservaciones
    fdp=zeros(1,nesp+1);
    x2=X(k,:);
    for iesp=1:nesp+1
        for in=1:ng
            gmd=gmdEsp{iesp};
            fdp(iesp)=fdp(iesp)+gmd.PComponents(in)*...
                mvnpdf(x2,gmd.mu(in,:),gmd.Sigma(:,:,in));
        end %de gaussiana
    verosimilitud(iesp)=log10(max(fdp(iesp),realmin));
    end %de especie
    [~,espopt]=max(verosimilitud);
    Z(k)=espopt;
    ProbZ(k,:)=verosimilitud;
end %de frame
kk=1;


%
% Función que obtiene el clasificador por distancia mínima
%
function clasificador=fitDistMin(X,Y)
% X: matriz de predictores (observación,parámetro)
% Y: vector de clases (observación)
% clasificador: clasificador; vector (clases)

global PATRON;
nesp=PATRON.nesp;
npar=size(X,2);
mu=zeros(nesp+1,npar);
sigma=zeros(nesp+1,npar);
for iesp=1:nesp+1
    Xesp=X(Y==iesp,:);
    mu(iesp,:)=median(Xesp,1);
    sigma(iesp,:)=iqr(Xesp,1);
end
clasificador={mu,sigma};


%
% Función que clasifica mediante Distancia Mínima
%
function [Z,ProbZ]=DistMinPredict(clasificador,X)
% X: matriz de predictores (observación,parámetro)
% clasificador: clasificador
% Z: vector de predicciones (observación)
% ProbZ: Probabilidad de las predicciones: matriz (observación,especie)

global PATRON;
nesp=PATRON.nesp;
nobservaciones=size(X,1);
mu=clasificador{1};
Z=zeros(nobservaciones,1);
ProbZ=zeros(nobservaciones,nesp+1);
distancia=zeros(nesp+1,1);
rangopar=max(X,[],1)-min(X,[],1);

for k=1:nobservaciones
    x2=X(k,:);
    for iesp=1:nesp+1
        mux=mu(iesp,:);
        distpar=abs(x2-mux);
        distnorm=distpar./rangopar;
        distancia(iesp)=sqrt(sum(distnorm.^2));
    end %de especie
    
    score=sum(distancia)./distancia;
    ProbZ(k,:)=score./sum(score);
    [~,espmin]=min(distancia);
    Z(k)=espmin;
end %de frame
kk=1;