% Hald example
clear all

load hald

Ngrid=100;
[N,D]=size(ingredients);
xmin=min(ingredients); xmax=max(ingredients);
xgrid=zeros(Ngrid,D);
for i=1:D
    xgrid(:,i)=linspace(xmin(i),xmax(i),Ngrid);
end


labs={'tricalcium aluminate','tricalcium silicate','tetracalcium aluminoferrite','beta-dicalcium silicate'};

figure(2)
set(gcf,'Name','Grafico')
clf;
leg={'Dati','Fit retta','Fit ordine 2','Fit 4 variabili'};
h=zeros(4,1);

% Fit a 4 variabili
    mat=[ones(N,1) ingredients];
    fprintf('Fit a 4 variabili: ')
    coeffb4=mat\heat;
    fprintf('Coefficienti: %f %f %f %f %f - ',coeffb4)
    fitb4=mat*coeffb4;
    errorb4=sqrt(sum((fitb4-heat).^2)/N);
    fprintf('Errore: %6.2f - ',errorb4)
    [r,p]=corrcoef(ingredients(:,i),heat);
    fprintf('Coeff. correl.: %f (p-value: %f)',r(1,2),p(1,2));
fprintf('\n')


for i=1:D
    fprintf('\n*** Variable %s\n',labs{i});
    subplot(2,2,i)
    h(1)=plot(ingredients(:,i),heat,'+b');
    hold on
    xlabel('Percentage')
    ylabel('Heat [cal/gm]')
    title(labs{i})
    
    % Fit minimi quadrati
    % ordine 1
    mat=[ones(N,1) ingredients(:,i)];
    fprintf('Fit con retta: ')
    coeff1=mat\heat;
    fprintf('Coefficienti: %f %f - ',coeff1)
    fit1=mat*coeff1;
    error1=sqrt(sum((fit1-heat).^2)/N);
    fprintf('Errore: %6.2f - ',error1)
    [r,p]=corrcoef(ingredients(:,i),heat);
    fprintf('Coeff. correl.: %f (p-value: %f)',r(1,2),p(1,2));
    fprintf('\n')
    fit1=[ones(Ngrid,1) xgrid(:,i)]*coeff1;
    h(2)=plot(xgrid(:,i),fit1,'r');
    
    % ordine 2
    mat=[ones(N,1) ingredients(:,i) ingredients(:,i).^2];
    fprintf('Fit del secondo ordine: ')
    coeff2=mat\heat;
    fprintf('Coefficienti: %f %f %f - ',coeff2)
    fit2=mat*coeff2;
    error2=sqrt(sum((fit2-heat).^2)/N);
    fprintf('Errore: %6.2f - ',error2)
    fprintf('\n')
    fit2=[ones(Ngrid,1) xgrid(:,i) xgrid(:,i).^2]*coeff2;
    h(3)=plot(xgrid(:,i),fit2,'g');
    
    [sorted_ingredients ind]=sort(ingredients(:,i));
    h(4)=plot(sorted_ingredients,fitb4(ind),'c');
    
    legend(h,leg,'Location','Best')
    
end