%% LABORATORIO 2 - METODI STATISTICI PER LA BIOINGEGNERIA
%% Canale B - A.A. 2024/2025
%% Docente: Enrico Longato

%% Parte 1
clc        
clear all  
close all  

%% Caricare i dati
load data_IVGTT % Osservare il workspace dopo il caricamento

%% Estrarre le variabili ed etichette di interesse
subset_labels = {'glucose_ss', 'insulin_basal'}; % Etichette in un cell array
% Per casa: provare a fare un vettore NON cell array
% ['glucose_ss', 'insulin_basal'] e rendersi conto di perché 
% l'operazione fallisce.

% Recupero gli indici
indices_of_subset_labels = []; % Inizializzo un array vuoto
for label_curr = subset_labels
    % Cerco la stringa
    i_label_curr = find(strcmp(label_curr, labels));
    
    indices_of_subset_labels = [indices_of_subset_labels, i_label_curr];
    
end

% Salvo come richiesto dalla consegna
subset_data = data(:, indices_of_subset_labels);

% Alternativamente (all'esame è consentito), potevo ispezionare il cell
% array manualmente e "hard codificare" la selezione.
% subset_labels = labels([7 11]);
% subset_data = data(:, [7 11]);

%% Aggiungere le unità di misura
subset_units = {'mg/dL', 'pmol/mL'}; % Stesso ordine
for i = 1:length(subset_labels)
    subset_labels{i} = [subset_labels{i} ' (' subset_units{i} ')'];
end

% Si poteva usare l'istruzione strcat invece del ciclo for
% subset_labels = ...
%       strcat(subset_labels, {' (', ' ('}, subset_units, {')', ')'})

%% Visualizzare uno scatter plot delle due variabili
% Apro una figura vuota (se voglio, posso assegnarle un ID con figure(ID))
figure
scatter(subset_data(:, 1), subset_data(:, 2)) % Scatter plot senza fronzoli
xlabel(subset_labels{1}) 
ylabel(subset_labels{2})

% Non abbiate paura che viene una "s" strana come se fosse un pedice;
% non è un problema per l'esame, in cui si valuta il contenuto e non certo
% l'estetica.

%% Visualizzare i due boxplot in due figure separate
for k = 1:size(subset_data, 2) % Numero di colonne di subset_data
    % Apro una figura per iterazione e, per non confondermi, assegno un ID
    % esplcito a ciascuna
    figure(10 + k) % 10 è un ID qualunque di partenza "che non dà fastidio"

    boxplot(subset_data(:, k))
    xlabel(subset_labels{k})

end

%% Visualizzare i due boxplot in due subplot
% Apro la figura fuori dal ciclo, perché è una sola
figure(20) % Altro ID "che non dà fastidio" perché "distante" da 10
for k = 1:size(subset_data, 2)
    % L'istruzione qui sotto si legge:
    % In un gruppo di subplot di 1 riga e 2 colonne, 
    % disegna nel k-esimo subplot
    subplot(1, 2, k) 

    boxplot(subset_data(:, k)) 
    xlabel(subset_labels{k})
    ylabel(subset_units{k})

end

%% Parte 2

%% Calcolare media, moda e mediana delle due variabili
stats_data = struct(); % Inizializzo la struct

for k = 1:size(subset_data, 2)
    % Calcoliamo le statistiche
    mean_curr = mean(subset_data(:, k));
    mode_curr = mode(subset_data(:, k)); % Chiedersi se ha senso per v. continue
    median_curr = median(subset_data(:, k));

    % Inseriamo nella struct in maniera intuitiva: MATLAB inizializza
    % "da solo" i campi come campi e i vettori come vettori
    stats_data.mean(k) = mean_curr; % MATLAB capisce cosa intendevamo
    stats_data.mode(k) = mode_curr; 
    stats_data.median(k) = median_curr; 
end

%% Parte 3

%% Uso del find 

% Trovare glucose_ss maggiore di 90 mg/dL con find
i_gss_high = find(subset_data(:, 1) > 90);

% Mostrare i boxplot della sottopopolazione
figure(30) % Solito ID "che non dà fastidio" perché "distante" dagli altri
for k = 1:size(subset_data, 2)
    subplot(1, 2, k) 

    boxplot(subset_data(i_gss_high, k))
    xlabel(subset_labels{k})
    ylabel(subset_units{k})

end

%% Uso diretto dell'array booleano (viene lo stesso plot)
i_gss_high_mask = subset_data(:, 1) > 90; % Sarebbe l'argomento del find
% Per casa: confrontare i_gss_high e i_gss_high_mask dal workspace

figure(31) 
for k = 1:size(subset_data, 2)
    subplot(1, 2, k) 

    boxplot(subset_data(i_gss_high_mask, k)) % Mask qui invece degli indici!
    xlabel(subset_labels{k})
    ylabel(subset_units{k})

end

%% Parte 4

%% Istogrammi con 20 bin
N_bins = 20;

for k = 1:size(subset_data, 2)
    % Si parte dall'istogramma delle frequenze assolute
    [counts, centres] = hist(subset_data(:, k), N_bins); 
    % NB: hist è sufficiente per i nostri scopi, ma all'esame potete usare
    %     quello che volete; basta che la figura sia giusta.
    
    figure(40)
    subplot(2, 2, k)
    bar(centres, counts)
    xlabel(subset_labels{k})
    ylabel('Frequenze assolute')

    % Per trovare le frequenze relative, dividiamo per il numero di
    % campioni
    freq_rel = counts / size(subset_data, 1);

    subplot(2, 2, k + 2)
    bar(centres, freq_rel)
    xlabel(subset_labels{k})
    ylabel('Frequenze relative')

end

%% Gaussianità "a occhio"
disp([subset_labels{1} ' potrebbe essere gaussiana.'])
disp([subset_labels{2} ' è troppo asimmetrica per essere gaussiana.'])