% 16 scrittori di romanzi sono stati valutati da un campione di lettori che
% hanno espresso opinioni sul tipo di contenuto e sul modo in cui sono state
% scritte le loro opere in base alle seguenti variabili 
% X1: giudizio complessivo, 
% X2: leggibilità, 
% X3: politica, 
% X4: fantasia, 
% X5: desiderio di rilettura, 
% X6: attualità. 
% 
% Importare i dati delle 6 variabili dal file di Excel denominato
% scrittori.xlsx utilizzando la funzione readtable
%
% QUESITI
% 1) Rappresentare i dati in termini di icone (o facce) 
% in modo da individuare eventuali scrittori anomali.
% 2) Calcolare la matrice di correlazione. Commentare il coefficiente di
% correlazione lineare tra le variabili X5 e X2. Costruire e commentare il
% diagramma di dispersione tra le variabili X5 (in asse delle ordinate) e
% X2 (in asse delle ascisse)  e commentare la relazione. Illustrare se in
% questo diagramma di dispersione c'è qualche valore anomalo bivariato.
% 3) Effettuare una riduzione delle dimensioni (operare sulle variabili in
% forma standardizzata) 
% 4)  Scegliere un numero di CP appropriato
% 5) Calcolare e commentare le correlazioni tra CP e le variabili
% originarie. Calcolare e interpretare le comunalità.
% 6) Interpretare le CP estratte. Che caratteristiche hanno gli scrittori
% che presentano valori elevati della prima CP? Che caratteristiche hanno
% gli scrittori che presentano valori elevati della seconda CP?
% Costruire il biplot inserendo come punti riga le componenti principali
% standardizzate Z*V*\Lambda^-0.5 e come punti colonna (frecce) le correlazioni.
% 7) Costruire il biplot tramite la svd utilizzando omega=1 e alpha =0,
% ossia punti riga come \sqrt(n-1)*U (componenti principali standardizzate)
% e punti riga come Gamma*V (correlazioni). Veficare che i punti riga e
% colonna calcolati nei due metodi 6) e 7) coincidono
%
% I quesiti che seguono possono essere risolti solo dopo aver svolto la
% cluster analysis
%
%
% 8) Applicare ai dati standardizzati la cluster analysis basata sul metodo
% delle k-medie utilizzando 3 gruppi. 
% 9) Calcolare i centroidi dei 3 gruppi in termini di variabili
% originarie e di variabili in forma standardizzata.  Proiettare i
% centroidi nello spazio delle prime due componenti principali.
% Aggiungere al biplot precedente i 3 centroidi. Interpretare
% i 3 gruppi che sono stati trovati.
% 10) Calcolare la distanza Euclidea e la distanza di Mahalanobis di ogni
% scrittore dal centroide. Calcolare e commentare il coefficiente di
% correlazione tra le due graduatorie.


%% Soluzione
% Importare i dati delle 6 variabili dal file di Excel denominato
% scrittori.xlsx utilizzando la funzione readtable
Xtable=readtable('scrittori.xlsx','ReadRowNames',1,'Range','A1:G17');


%% 1) Rappresentare i dati in termini di icone (o facce) 
% in modo da individuare eventuali scrittori anomali.
RowNam=Xtable.Properties.RowNames;
ColNam=Xtable.Properties.VariableNames;
X=table2array(Xtable);
% rappresentazione tramite stelle
figure
glyphplot(X, 'glyph','star','obslabels',RowNam)
title('Rappresentazione tramite stelle')
figure
% rappresentazione tramite facce
glyphplot(X, 'glyph','face','obslabels',RowNam)
title('Rappresentazione tramite facce')

% Commento: Tolkien Rowling hanno una faccia 
% (stella) molto simile
% Lucarelli sembra avere un profilo simile 
% a Allende e King

%% 2) Calcolare la matrice di correlazione. Commentare il coefficiente di
% correlazione lineare tra le variabili X5 e X2. Costruire e commentare il
% diagramma di dispersione tra le variabili X5 (in asse delle ordinate) e
% X2 (in asse delle ascisse)  e commentare la relazione. Illustrare se in
% questo diagramma di dispersione c'è qualche valore anomalo bivariato.
correl=corr(X);
array2table(correl,'RowNames',ColNam,...
    'VariableNames',ColNam)

i=2; j=5;
plot(X(:,i),X(:,j),'o')
text(X(:,i),X(:,j),RowNam)
xlabel(ColNam(i),'Interpreter','none')
ylabel(ColNam(j),'Interpreter','none')

% Commento:
% Nasar sembra essere un valore anomalo bivariato. E' uno scrittore che
% presenta una leggibilità intorno alla media ma un desiderio di rilettura
% molto più basso degli altri


%% 3) Effettuare una riduzione delle dimensioni (operare sulle variabili in
% forma standardizzata)
[n,p]=size(X);
[Z,mu,Sigma]=zscore(X);

% Per trovare le componenti principali di seguito si utilizza la funzione di 
% MATLAB pca. 
[Coeff,Score,latent] =pca(Z);
% coeff sono gli autovettori (prima colonna = primo autovettore associato
% all'autovalore più grande ...)
% score è la matrice nx p che contiene le componenti principali (Z*V)
% latent è il vettore che contiene gli autovalori ordinati (dal più grande
% al più piccolo)
disp(cumsum(latent)/p)

% Le prime due componenti principali spiegano il 75.18% della varianza
% totale
% Adottando il criterio di considerare solo gli autovalori >1 dobbiamo estrarre le prime due componenti principali
% Adottando il criterio 0.95^p=0.69 posso limitarmi a considerare solo le
% prime due PC.
% La prima componente spiega il 46% circa della variabilitá totale
% presente nei dati, mentre il primo piano fattoriale (dato
% dalle prime due componenti principali) spiega circa il
% 75% della variabilitá.

% 5) Calcolare e commentare le correlazioni tra CP e le variabili
% originarie. Calcolare e interpretare le comunalità.

% MatrComp contiene le correlazioni tra le variabili e le prime due
% componenti principali
MatrComp=Coeff(:,1:2).*sqrt(latent(1:2)');


% 6) Interpretare le CP estratte. Che caratteristiche hanno gli scrittori
% che presentano valori elevati della prima CP? Che caratteristiche hanno
% gli scrittori che presentano valori elevati della seconda CP?
% Analisi delle correlazioni tra le variabili originarie e la prime 
% CP

% Se si adotta il segno che segue per la prima colonna di MatrComp
disp(array2table(-MatrComp(:,1),'RowNames',ColNam,'VariableNames',{'corr_withPC1'}))
% Osservazione: se si cambia di segno la prima colonna di MatrComp1 bisogna
% ricordarsi di cambiare di segno anche score(:,1)

% Gli scrittori che presentano valori elevati della prima CP  sono
% caratterizzati dai temi fantasy elevata leggibilità, temi non di
% attualità e bassi temi di politica

disp(array2table(MatrComp(:,2),'RowNames',ColNam,'VariableNames',{'corr_withPC2'}))
% Gli scrittori che presentano valori elevati della seconda CP hanno avuto
% un alto giudizio complessivo, elevata leggibilità e  ispirano alla
% rilettura

% Se si adotta il segno che segue per la seconda colonna di MatrComp
disp(array2table(MatrComp(:,1:2),'RowNames',ColNam,'VariableNames',{'corr_withPC1', 'corr_withPC2'}))

% Comunalità
disp(array2table(sum(MatrComp.^2,2),'RowNames',ColNam,'VariableNames',{'Comun'}))
% La variabile spiegata meglio è l’attualità

% Rappresentare i punti riga come le componenti principali standardizzate
% Cambio segno a MatrComp(:,1) e di conseguenza anche a 
MatrComp(:,1)=-MatrComp(:,1);
Score(:,1)=-Score(:,1);
% Punti riga rappresentati dalle componenti principali in forma
% standardizzata
PuntiRiga=zscore(Score);
figure
plot(PuntiRiga(:,1),PuntiRiga(:,2),'o')
text(PuntiRiga(:,1),PuntiRiga(:,2),Xtable.Properties.RowNames)

hold('on')
PuntiColonna=MatrComp;
zeroes = zeros(p,1); 
quiver(zeroes,zeroes,PuntiColonna(:,1),PuntiColonna(:,2))
varlabs=Xtable.Properties.VariableNames;
dx=0;
dy=0;
text(PuntiColonna(:,1)+dx,PuntiColonna(:,2)+dy,varlabs,'Color','b','Interpreter','none');

axislim=axis;
% Aggiungo l'asse x
line([axislim(1);axislim(2)], [0;0], 'Color','black');
% Aggiungo l'asse y
line([0;0],[axislim(3);axislim(4)], 'Color','black');
title('Biplot costruito senza utilizzare la funzione svd')


%% Biplot (passando attraverso la svd)
[U,Gammastar,V]=svd(Z,'econ');
sqn1=sqrt(n-1);
Gamma=Gammastar/sqn1;
omega=1;
alpha=0;
PuntiRigaSVD=sqn1^omega*U(:,1:2)*Gamma(1:2,1:2)^alpha;
PuntiColonnaSVD=V(:,1:2)*(Gamma(1:2,1:2)^(1-alpha))*sqn1^(1-omega);

plot(PuntiRigaSVD(:,1),PuntiRigaSVD(:,2),'o')
text(PuntiRigaSVD(:,1),PuntiRigaSVD(:,2),Xtable.Properties.RowNames)

zeroes = zeros(p,1); 
quiver(zeroes,zeroes,PuntiColonnaSVD(:,1),PuntiColonnaSVD(:,2))
varlabs=Xtable.Properties.VariableNames;
dx=0.02;
dy=0.03;
text(PuntiColonnaSVD(:,1)+dx,PuntiColonnaSVD(:,2)+dy,varlabs,'Color','b','Interpreter','none');

axislim=axis;
% Aggiungo l'asse x
line([axislim(1);axislim(2)], [0;0], 'Color','black');
% Aggiungo l'asse y
line([0;0],[axislim(3);axislim(4)], 'Color','black');

% Commenti:
% Angolo piccolo tra politic e attual==> le due  variabili sono fortemente
% correlate ==>es. leggibilità e rilettura
%......................
title('Biplot costruito utilizzando la funzione svd')

% Verifico che i punti riga (componenti principali standardizzate
% calcolati come Z*V*\Lambda^(-1/2) (matrice PuntiRiga) siano uguali ai
% quelli calcolati tramite la svd (matrice PuntiRigaSVD)
disp('Controllo sui punti riga trovati con i due metodi')
max(abs(PuntiRiga(:,1:2)-PuntiRigaSVD),[],'all')

% Verifico che i punti colonna (correlazioni
% V*\Lambda^(1/2) (matrice PuntiColonna) siano uguali ai
% quelli calcolati tramite la svd (matrice PuntiColonnaSVD)
max(abs(PuntiColonna(:,1:2)-PuntiColonnaSVD),[],'all')


%% La parte di seguito è riferita alla cluster analysis

% 8) Applicare ai dati standardizzati la cluster analysis basata sul metodo
% delle k-medie utilizzando 3 gruppi. 
[idx,C]=kmeans(Z,3,'Replicates',100);
% spmplot(X,idx)


% 9) Calcolare i centroidi dei 3 gruppi in termini di variabili
% originarie e di variabili in forma standardizzata.  Proiettare i
% centroidi nello spazio delle prime due componenti principali.
% Aggiungere al biplot precedente i 3 centroidi. Interpretare
% i 3 gruppi che sono stati trovati.

disp('Centroidi in forma di scostamenti standardizzata')
disp(C)
disp('Centroidi in termini dei dati originari')
disp(C.*Sigma+mu)

% Controllo di aver fatto i calcoli giusti
% Centroide del gruppo 1 
mean(X(idx==1,:))

% Proiezione dei centroidi nello spazio delle componenti principali
Centrproiettati=C*V(:,1:2)*(sqn1^(omega-1))*Gamma(1:2,1:2)^(alpha-1);
plot(Centrproiettati(:,1),Centrproiettati(:,2),'k+','MarkerSize',12,'LineWidth',2)
sc=size(C,1);
centrlabels=[repmat('C',sc,1) num2str((1:sc)')];
text(Centrproiettati(:,1),Centrproiettati(:,2)+0.04,centrlabels);

% text(Centrproiettati(:,1),Centrproiettati(:,2),'+','FontSize',14)
% Aggiunta al grafico del gruppo di appartenenza
text(PuntiRiga(:,1)-0.03,PuntiRiga(:,2),num2str(idx))

%% 10) Calcolare la distanza Euclidea e la distanza di Mahalanobis di ogni
% scrittore dal centroide. Calcolare e commentare il coefficiente di
% correlazione tra le due graduatorie.
mu=mean(X);
% Calcolo delle distanze euclidee e di Mahalanobis dal centroide generale
dMahalanobis=zeros(n,1);
dEuclidea=dMahalanobis;
invSigma=inv(cov(X));
for i=1:n
dMahalanobis(i) = sqrt((X(i,:)-mu)*invSigma*(X(i,:)-mu)');
dEuclidea(i)= sqrt((X(i,:)-mu)*(X(i,:)-mu)');
end
figure
corr(dMahalanobis,dEuclidea,'type','Spearman')
plot(dMahalanobis,dEuclidea,'o')
% Commento la correlazione tra i posti d'ordine è scarsa.