%mreg.m

%   Copyright 2017 Rafael Almar (IRD, France)- rafael.almar@ird.fr
function [ beta, varbeta, sigma2, residuals, preds ,R2] = mreg(y, X)
% MREG  Performs a multiple regression fit given y and X
% USE:  [ beta, varbeta, sigma2, residuals, preds ] = mreg(y, X)

%Attention; la première colonne doit ^tre des 1 pr calculer Bo
%beta= linear coefficeint Y=B0 + B1X1 + B2X2 ...
%sigma2= variance of the error term (sigma^2)
%R2 :It is the proportion of variability in a data set that is accounted
%for by the statistical model. It provides a measure of how well future
%outcomes are likely to be predicted by the model.

%CHECK INPUT ARGUMENTS
if ( nargin ~= 2 )
   error('Need two arguments: y and X.');
end

%ERROR CHECKS
detX = det( X' * X );

if ( abs(detX) < 0.000000001 )
   disp(['Determinant = ',num2str(detX)]);
   error('X is singular, or almost singular.');
end
%END ERROR CHECKS


if ~isempty(X),
   beta = inv( X' * X ) * ( X' * y );
   %or:  beta = ( X' * X )^(-1) * ( X' * y );
   preds = X * beta;
   residuals = y - preds;
else
   beta = 0;
   preds = zeros(size(y));
   residuals = y;
end;

%An unbiased estimate of sigma^2
sigma2 = sum( residuals.^2 ) / ( length(y) - length(beta) );

varbeta = sigma2 * inv( X' * X );

%Calcul des coefficients de détermination

for i=1:length(X(1,:))
   R =( corrcoef(y,(beta(i).*X(:,i)'))).^2;
    R2(i)=R(2,1);
end