%exercice 3
close all force;clear all;

disp('Ouverture du fichier de données...');
load data_exercice3
%Essai pour un sujet
% [p,tbl,stats1,terms1]=anovan(Y(:,2),{g1 g2 g3 g4 g5 g6 g7 g8},'model','linear','sstype',3,'varnames',varnames);
% %coeffficients du modèle
% for i=1:length(stats1.coeffs)
%     disp([cell2mat(stats1.coeffnames(i)),' :  ',num2str(stats1.coeffs(i))]);
% end
% return
utilites=zeros(25,30);
importances=zeros(8,30);
%indices inférieur et supérieur des niveaux correspondant à chaque facteur
%dans les grands vecteurs de paramètre (1 = constante du modèle);
numinf=[2 5 8 11 13 16 19 23];
numsup=[4 7 10 12 15 18 22 25];
disp('Calcul des modèles individuels...');
for n=1:30
    [p,tbl,stats1,terms1]=anovan(Y(:,n),{g1 g2 g3 g4 g5 g6 g7 g8},'model','linear','sstype',3,'varnames',varnames);
    Rsquared(n)=1-cell2mat(tbl(size(tbl,1)-1,2))/cell2mat(tbl(size(tbl,1),2));
    utilites(:,n)=stats1.coeffs;
    for i=1:8
        amplitude(i)=max(utilites(numinf(i):numsup(i),n))-min(utilites(numinf(i):numsup(i),n));
    end
    for i=1:8
        importances(i,n)=amplitude(i)/sum(amplitude);
    end
    close all force
end
%graphe des coefficients de détermination
figure;plot(Rsquared,'k*-');grid on;xlabel('participant');ylabel('coeff. determ.');
%graphe des utilités moyennes
utilites_moyennes=mean(utilites,2);
figure;
for i=1:8
    plot([numinf(i):numsup(i)],utilites_moyennes(numinf(i):numsup(i)),'k*-');hold on;
end

%boxplot des importances
figure;
boxplot(importances',varnames);ylabel('importances');xlabel('facteurs');

%on enlève les sujets pour lesquels R2<0.7
limite=0.7;
sujets_ok=find(Rsquared>=limite);
sujets_nonok=find(Rsquared<limite);
disp(['Limite : ',num2str(limite)]);
s='on enlève les participants n°';
for i=1:length(sujets_nonok)
    s=[s,int2str(sujets_nonok(i)),' '];
end
disp(s);
importances_ok=importances(:,sujets_ok(:));
utilites_ok=utilites(:,sujets_ok(:));
%evaluations individuelles du produit A et du produit B
modA=[1 2 5 8 11 13 16 19 23];
evalA=sum(utilites_ok(modA(:),:));
modB=modA;modB(2:length(modB))=modB(2:length(modB))+1;
evalB=sum(utilites_ok(modB(:),:));
prefA=100*length(find(evalA>evalB))/28;

%classification des participants
disp('Classification des participants...');
Z=linkage((importances_ok)','ward');
figure;dendrogram(Z);
%coupure de l'arbre
nbre_classes=2;
c = cluster(Z,'maxclust',nbre_classes);


figure;
for i=1:nbre_classes
    subplot(nbre_classes,1,i);toto=find(c==i);
    z=importances_ok(:,toto(:));
    boxplot(z',varnames);ylabel('importances');xlabel('facteurs');
    title(['groupe ',int2str(i),' (',int2str(length(toto)),' pers.)']);
end
figure;
symbole={'k*-','ro-','bd-'};
for n=1:nbre_classes
    toto=find(c==n);
    for i=1:8
        plot([numinf(i):numsup(i)],mean(utilites_ok(numinf(i):numsup(i),toto(:)),2),char(symbole(n)));hold on;
    end
end