[nz,ntv] = size(obs);
% means
disp(' ')
mobs = mean(obs);
disp (['mean observed values = ' num2str(mobs)])
mest = mean(est);
disp (['mean estimated values = ' num2str(mest)])
% standart deviations
sdobs = std(obs);
disp (['std dev observed values = ' num2str(sdobs)])
sdest = std(est);
disp (['std dev estimated values = ' num2str(sdest)])
% assessment of goodness-of-fit
% correlation and coefficient of determination (explained variance)
[corre,p] = corrcoef([obs est]);R2=corre.^2;
disp (['R2 corre^2 = ' num2str(R2(1,2)) ' p = ' num2str(p(1,2))])
% mean absolute error
imae = sum(abs(obs-est))/ntv;
disp (['mae = ' num2str(imae)])
% root mean square error
rmse = sqrt(sum((obs-est).^2)/ntv);
disp (['rmse = ' num2str(rmse)])
disp(' ')
% Nash-Sutcliffe coeff of efficiency
E = 1-(sum((obs-est).^2)/sum((obs-mobs).^2));
disp (['N-Sut efficiency (E) = ' num2str(E)])
stat = bootstrp(500,@Nash_Sut,obs,est);
[a b]=ci95(stat);
disp ([ 'E 95% ci = ' num2str(a) ' ' num2str(b)])
disp (' ');clear stat
% log Nash-Sutcliffe coeff of efficiency
Elog = 1-(sum((log(obs)-log(est)).^2)/sum((log(obs)-log(mobs)).^2));
disp (['Elog = ' num2str(Elog)])
stat = bootstrp(500,@log_Nash_Sut,obs,est);
[a b]=ci95(stat);
disp ([ 'E 95% ci = ' num2str(a) ' ' num2str(b)])
disp (' ');clear stat
% Modified Nash-Sutcliffe coeff of efficiency
E1 = 1-(sum(abs(obs-est))/sum(abs(obs-mobs)));
disp (['N-Sut efficiency modified (E1) = ' num2str(E1)])
stat = bootstrp(500,@m_Nash_Sut,obs,est);
[a b]=ci95(stat);
disp ([ 'E1 95% ci = ' num2str(a) ' ' num2str(b)])
disp (' '); clear stat
% index of agreement
d = 1-(sum((obs-est).^2)/sum((abs(obs-mobs)+abs(obs-mobs)).^2));
disp (['Index of agreement (d) = ' num2str(d)])
stat = bootstrp(500,@I_agree,obs,est);
[a b]=ci95(stat);
disp ([ 'd 95% ci = ' num2str(a) ' ' num2str(b)])
disp (' '); clear stat
% Modified index of agreement
d1 = 1-(sum(abs(obs-est))/sum((abs(obs-mobs)+abs(obs-mobs))));
disp (['Index of agreement modified (d1) = ' num2str(d1)])
stat = bootstrp(500,@m_I_agree,obs,est);
[a b]=ci95(stat);
disp ([ 'd1 95% ci = ' num2str(a) ' ' num2str(b)])
disp (' '); clear stat