% ChiSquared Demo

% Brian D. Jeffs
% 10/19/05

echo on
X = randn(100000,1);

plot(X(1:1000),'+');
title('1000 Gaussian samples');
pause

Nbins = 100
Nsamps = 100000
[pdf_hat, x] = hist(X(1:Nsamps),Nbins);
binWidth = x(3)-x(2);
pdf_hat = pdf_hat/Nsamps/binWidth;
xSamps = [-5:.1:5];
trueGaussian = normpdf(xSamps);
plot(x,pdf_hat,'.',xSamps,trueGaussian,'r');
legend(['Histogram for ',num2str(Nsamps),' samples'],'True Gaussian pdf');
hold on
bar(x,pdf_hat);
title('Comparison of histogram pdf estimate to true pdf');
hold off
pause

Chi = X.^2;
nu = 1;
plot(Chi(1:1000),'+');
title('1000 /chi^2 samples, 1 degree of freedom');
pause

Nbins = 100
Nsamps = 100000
[pdf_hat, x] = hist(Chi(1:Nsamps),Nbins);
binWidth = x(3)-x(2);
pdf_hat = pdf_hat/Nsamps/binWidth;
xSamps = [-5:.1:5];
trueChi2 = chi2pdf(xSamps,nu);
plot(x,pdf_hat,'.',xSamps,trueChi2,'r');
legend(['Histogram for ',num2str(Nsamps),' samples'],['True \chi^2 pdf, ',num2str(nu),' degrees of freedom']);
hold on
bar(x,pdf_hat);
title('Comparison of histogram pdf estimate to true pdf');
hold off
pause

X = randn(100000,1);
Chi = Chi + X.^2;  % Adding two squared Gaussians
nu = 2;

Nbins = 100
Nsamps = 100000
[pdf_hat, x] = hist(Chi(1:Nsamps),Nbins);
binWidth = x(3)-x(2);
pdf_hat = pdf_hat/Nsamps/binWidth;
xSamps = [-5:.1:5];
trueChi2 = chi2pdf(xSamps,nu);
plot(x,pdf_hat,'.',xSamps,trueChi2,'r');
legend(['Histogram for ',num2str(Nsamps),' samples'],['True \chi^2 pdf, ',num2str(nu),' degrees of freedom']);
hold on
bar(x,pdf_hat);
title('Comparison of histogram pdf estimate to true pdf');
hold off
pause

echo off