Thread: Help with Matlab - advanced stats

Hi,

I need some help with a matlab m-file I was given to help me do some statistical calculations. I am not too familiar with Matlab and am having problems adapting the m-file below to do what I need it to.

The specific help I need is:

1. adapt the m-file to take 30x23 matrices (.csv files)
2. debug, I kept getting error messages but that could be purely due to my inept coding
3. and last but not least how to run it, or is it as simple as pressing f5?

I would appriciate any help on this. Here is the code as given to me:

function [mobs,mest,sdobs,sdest,R2,imae,rmse,E,Elog,E1,d,d1] = goodfit2(obs,est)

% this is a function I developed to calculate
% different indices for checking the goodness-of-fit
% of any model (not valid for binary outputs)
%
% inputs are:
% est - estimated values by the model (should be a vector (ntv,1))
% obs - values the model should be compared to (should be a vector (ntv,1))
%
% outputs are:
% R2 - square of the correlation btw observed and estimated values
% MAE - mean absolut error
% RMSE - root mean square error
% E - coefficient of efficiency - Nash and Sutcliffe (1970)
% Elog - log Nash-Sutcliffe coefficient of efficiency
% E1 - modified Nash-Sutcliffe coefficient of efficiency
% d - index of agreement (Willmott 1981)
% d1 - modified index of agreement


[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

Originally Posted by Kean

Hi,

I need some help with a matlab m-file I was given to help me do some statistical calculations. I am not too familiar with Matlab and am having problems adapting the m-file below to do what I need it to.

The specific help I need is:

1. adapt the m-file to take 30x23 matrices (.csv files)
2. debug, I kept getting error messages but that could be purely due to my inept coding
3. and last but not least how to run it, or is it as simple as pressing f5?

I would appriciate any help on this. Here is the code as given to me:

What error messages?

It does not need adaptation it does not care what the dimensions of the data matrix are.

Save the file to a directory on the search path and call it goodfit2.m

Import the data using the data import tool. The use the command:

[mobs,mest,sdobs,sdest,R2,imae,rmse,E,Elog,E1,d,d1] = goodfit2(INobs,INest)

where INobs are the imported observations, and INest are the imported estimates.

.

Thanks Constantine, I forgot to mention that each column is a separate data set ie. hydrology station. Will this m-file return values for each column/set separately?

??? Error using ==> minus
Matrix dimensions must agree.

Error in ==> goodfit2 at 54
E = 1-(sum((obs-est).^2)/sum((obs-mobs).^2));