[Fat Age Wt Ht Neck Chest Waist Hip Thigh Knee Ankle Bicep Forearm Wrist] ...
    = textread('bodyfat.txt','%f%d%f%f%f%f%f%f%f%f%f%f%f%f','headerlines',1);

x = Waist;
y = Fat;
n = length(y);
b1 = corr(x,y)*std(y)/std(x)
b0 = mean(y) - b1*mean(x)
'Summary'
quantile(y-(b0+b1*x),[0, 0.25, 0.5, 0.75, 1])

figure(1);
plot(x,y,'k+');
hold on;
plot([min(x) max(x)], [b0+b1*min(x) b0+b1*max(x)], 'b-');
hold off;

figure(2);
plot(x, y-(b0+b1*x), 'ko');

% Bootstrap...

for j=1:10000
    idx = randi(n,n,1);
    x = Waist(idx);
    y = Fat(idx);
    b1 = corr(x,y)*std(y)/std(x);
    b0 = mean(y) - b1*mean(x);
    b0_list(j) = b0;
    b1_list(j) = b1;
end

b0_bar = mean(b0_list)
b0_std = std(b0_list)
b1_bar = mean(b1_list)
b1_std = std(b1_list)

figure(2); hist(b0_list, 300);
figure(3); hist(b1_list, 300);
quantile(b0_list, [0.025, 0.975])
quantile(b1_list, [0.025, 0.975])