function [a,R,N] = cronbach_stand(X)
% CRONBACH_STAND Standardized Cronbach's Alpha
%   [a,R,N]=cronbach_stand(X) calculates a standardized (respectively
%   ranked; note the text below) Cronbach's alpha of the data set X.
%
%   a is the standardized Cronbach's alpha.
%   R is the upper triangular Spearman inter-correlation matrix
%   N is the number of valid items
%   X is the data set. Columns are items.
%
%   If there are items with zero variance, these items are excluded from
%   the calculation. There is a warning in this case; alpha might be
%   guessed to high.
%
%   Frederik Nagel, 2006, frederik.nagel@hmt-hannover.de
%
%
%   This is a corrected version of the function "cronbach.m", written by 
%   Frederik Nagel in 2006. There was one major bug in "cronbach.m" in line 
%   69, eventually resulting in a too small number of items and therefore 
%   in a too small alpha-coefficient. A minor bug was in line 60, resulting 
%   in a somewhat misleading info-output.
%
%   Standardized Cronbach's Alpha is a derivate of the original alpha, 
%   presented by Cronbach 1951 (Cronbach, Lee J. (1951): Coefficient alpha 
%   and the internal structure of tests. Psychometrika 16(3), p. 297-334).
%   The there mentioned alpha-formula is implemented e.g. by my function 
%   "cronbach_alpha.m".
%
%   Whereas Cronbach's alpha-formua (and therefore the function 
%   "cronbach_alpha.m") is not able to cope with input data (X)
%   that include missing values (except by excluding the whole object 
%   providing one or more missing values), "cronbach_stand.m" works as well
%   with this kind of data (taking into account all available data).
%
%   cronbach_stand is equivalent to an "original"-alpha-coefficient 
%   computation (e.g. "cronbach_alpha.m") on the column-wise z-transformed 
%   (or "standardized") data X.
%   If the option "pearson" is changed to "spearman" this function yields
%   the same value as the "original"-alpha-coefficient computation (e.g.
%   "chronbach_alpha.m") on the column-wise rank-scaled (tied ranks) data X.
%   This optional value of "corr.m" can also be set to "kendall". This
%   produces a derivate of Cronbach's Alpha that relies on the rank
%   correlation by Kendall.
%
%
%   Markus Brückl, TU Berlin, makUs@brYkl.de, September 22, 2008


if nargin<1 || isempty(X)
   error('You shoud provide a data set.');
else
   % X must be a 2-dimensional matrix with at least 2 columns (and 2 rows)
   if size(X,2)<2
      error('Invalid data set.');
   end
end
% Items
N=size(X,2);

% Entries of the upper triangular matrix
e=(N*(N-1)/2);

% Pearsons's correlation coefficient
R = corr(X,'rows','pairwise','type','pearson');

% Coefficients from upper triangular matrix
Rh = triu(R,1);

% Mean of correlation coefficients
r = sum(sum(Rh))/e;

% If there are columns with zero variance, these have to be excluded.
if(isnan(r))
    disp('There are columns with zero variance!');
    disp('These columns have been excluded from the calculation of alpha!');        
    disp([num2str(sum(sum(isnan(Rh)))) ' coefficients of ' num2str(((N*N)-N)/2) ' have been excluded.']);    

    % Correct number of coefficients to be averaged 
    e = e-sum(sum(isnan(Rh)));
    
    % corrected mean of correlation coefficients 
    r = nansum(nansum(Rh))/e;
    
    % corrected number of items
    N = N - sum(isnan(diag(R)));
end
% Formular for standardized alpha
a=(N*r)/(1+(N-1)*r);
