%% Testo
% 1) Mostrare la funzione di densità della v.c \chi^2_g  con g che varia da 1 a 20 nell'intervallo x [0 30]. Inserire il grafico per ciascun g in un pannello separato. Inserire nel titolo di ogni grafico l'indicazione dei rispettivi gradi di libertà. Impostare il minimo e il massimo dell'asse x pari a [0 30]
% 2) Calcolare in una v.c. \chi^2_g con g=12
% a) la probabilità di ottenere un valore inferiore a 11.3
% b) la probabilità di ottenere un valore compreso nell'intervallo [10 11]
% c) la probabilità di ottenere un valore superiore a 20
% d) il quantile 0.9 (ossia il valore che lascia alla sua destra una probabilità del 10%)
% 3) Estrarre 1000000 di numeri casuali dalla \chi^2_g con g=12.
% Calcolare le frequenze relative empiriche
% a) di ottenere un valore inferiore a 11.3
% b) di ottenere un valore compreso nell'intervallo [10 11]
% c) di ottenre un valore superiore a 20
% d) calcolare il quantile 0.9 (ossia il valore che lascia alla sua destra 10000 casi)
% con le probabilità ottenute al punto sopra
%% Soluzione
%% 1) Mostrare la funzione di densità della v.c \chi^2_g  con g che varia da 1 a 20 nell'intervallo x [0 30]. Inserire il grafico per ciascun g in un pannello separato. Inserire nel titolo di ogni grafico l'indicazione dei rispettivi gradi di libertà. Impostare il minimo e il massimo dell'asse x pari a [0 30]

% Studio della densità della v.c. Chi2 al variare dei gradi di libertà

x=0:0.1:30;
for j=1:20
    
    subplot(4,5,j)

    denschi2=chi2pdf(x,j);
    
    plot(x,denschi2)
    title(['$\chi^2_{' num2str(j) '}$'],'Interpreter','latex')
    xlim([0 30])
end

%% 2) Calcolare in una v.c. \chi^2_g con g=12
%% a) la probabilità di ottenere un valore inferiore a 11.3
%% b) la probabilità di ottenere un valore compreso nell'intervallo [10 11]
%% c) la probabilità di ottenere un valore superiore a 20
%% d) il quantile 0.9 (ossia il valore che lascia alla sua destra una probabilità del 10%)
%

g=12;
% Prob X^2<b
b=11.3;
disp(['Prob ( \chi^2_g(12) <' num2str(b) ' = ' num2str(chi2cdf(b,g)) ')'])
% Prob (a< X^2<b)
a=10;
b=11;
disp(['(' num2str(a)  '< Prob \chi^2_g(12) <' num2str(b) ') = ' num2str(chi2cdf(b,g)-chi2cdf(a,g))])

% Prob X^2>b
b=20;
disp(['Prob ( \chi^2_g(12) >' num2str(b) ' = ' num2str(1-chi2cdf(b,g)) ')'])

% quantile 0.9
disp(['Quantile 0.9 in una \chi^2_g(12) =' num2str(chi2inv(0.9,g))])

%% 3) Estrarre 1000000 di numeri casuali dalla \chi^2_g con g=12.
%% Calcolare le frequenze relative empiriche
%% a) di ottenere un valore inferiore a 11.3
%% b) di ottenere un valore compreso nell'intervallo [10 11]
%% c) di ottenre un valore superiore a 20
%% d) calcolare il quantile 0.9 (ossia il valore che lascia alla sua destra 100000 casi)
%% con le probabilità ottenute al punto sopra

% Genero numeri casuali dalla distribuzione chi2 con
% Osservazione: i risultati che seguono non sono riproducibili a meno di fissare un seed dei umeri casuali

% Per riprodurre i risultati occorre decommentre la riga che segue
% Fisso uno starting point casuale nel generatore di numeri casuali
% rng(1000)

n=1000000;
y=chi2rnd(12,n,1);

% Frequenza relativa empirica di valori inferiori a 11.3
% MODO 1: estrarre i numeri minori di 11.3 e vedere la lunghezza del vettore
freqrequested=length(y(y<11.3))/n;
disp(['Frequenza relativa empirica di valori inferiori a 11.3=' num2str(freqrequested)])
% MODO 2: fare direttamente la somma dei valori TRUE nel vettore Booleano y<11.3
freqrequestedCHK=sum(y<11.3)/n;

% controllo che il modo 1 sia esattamente uguale al modo 2
disp('Controllo sui modi di estrarre le frequenze')
isequal(freqrequested,freqrequestedCHK)

% Frequenza relativa empirica di valori compresi tra 10 e 11
freqrequested=length(y(y>10 & y<11))/n;
disp(['Frequenza relativa empirica di valori compresi tra 10 e 11=' num2str(freqrequested)])

% Frequenza relativa empirica di valori superiori a 20
freqrequested=length(y(y>20))/n;
disp(['Frequenza relativa empirica di valori superiori a 20=' num2str(freqrequested)])

% Calcolo del quantile empirico 0.9
quanrequested=prctile(y,90);
disp(['Quantile emprico 0.9 =' num2str(quanrequested)])