ConcentrationTest

ConcentrationTest - Test homogeneity of concentration parameters.

 Test whether the concentration parameters for N groups of circular data are equal.
 The data are assumed to have Von Mises distributions. If the median concentration
 is inferior to 1, a randomization procedure is used.
 See "Statistical Analysis of Circular Data" (Fisher, p. 131).

  USAGE

    [h,p] = ConcentrationTest(angles,group,alpha,nRandomizations)

    angles         angles in radians
    group          group number for each sample
    alpha          optional significance level (default = 0.05)
    n              optional number of randomizations (default = 1000)

  SEE

    See also Concentration, CircularMean, CircularVariance, CircularConfidenceIntervals.

0001 function [h,p] = ConcentrationTest(angles,group,alpha,nRandomizations)
0002 
0003 %ConcentrationTest - Test homogeneity of concentration parameters.
0004 %
0005 % Test whether the concentration parameters for N groups of circular data are equal.
0006 % The data are assumed to have Von Mises distributions. If the median concentration
0007 % is inferior to 1, a randomization procedure is used.
0008 % See "Statistical Analysis of Circular Data" (Fisher, p. 131).
0009 %
0010 %  USAGE
0011 %
0012 %    [h,p] = ConcentrationTest(angles,group,alpha,nRandomizations)
0013 %
0014 %    angles         angles in radians
0015 %    group          group number for each sample
0016 %    alpha          optional significance level (default = 0.05)
0017 %    n              optional number of randomizations (default = 1000)
0018 %
0019 %  SEE
0020 %
0021 %    See also Concentration, CircularMean, CircularVariance, CircularConfidenceIntervals.
0022 
0023 % Copyright (C) 2004-2017 by Michaël Zugaro
0024 %
0025 % This program is free software; you can redistribute it and/or modify
0026 % it under the terms of the GNU General Public License as published by
0027 % the Free Software Foundation; either version 3 of the License, or
0028 % (at your option) any later version.
0029 
0030 global n N r;
0031 
0032 if nargin < 3,
0033     alpha = 0.05;
0034 end
0035 if nargin < 4,
0036     nRandomizations = 1000;
0037 end
0038 isradians(angles);
0039 
0040 % Number of groups (r), number of samples in each group (n), total number of samples (N)
0041 r = max(group);
0042 N = length(angles);
0043 n = zeros(1,max(group));
0044 for i = 1:r,
0045     n(i) = sum(group==i);
0046 end
0047 if min(n) < 10,
0048     error('This test requires a minimum of 10 samples per group.');
0049 end
0050 
0051 % Median concentration parameter (kappa)
0052 for i = 1:r,
0053     k(i) = Concentration(angles(group==i));
0054 end
0055 kappa = median(k);
0056 
0057 % Circular mean for each group (mu)
0058 for i = 1:r,
0059     mu(i) = CircularMean(angles(group==i));
0060 end
0061 
0062 if kappa >= 1,
0063     % Statistics for the data
0064     fr = ComputeStatistics(angles,group,mu);
0065     % p level
0066     p = 1 - fcdf(fr,r-1,N-r);
0067     p = 2 * min(p, 1-p);
0068 else
0069     % Center all samples
0070     for i = 1:r,
0071         angles(group==i) = angles(group==i)-mu(i);
0072     end
0073     % Compute the statistics for the data
0074     fr = ComputeStatistics(angles,group,mu);
0075     % Compute the statistics for each randomization
0076     for i = 1:nRandomizations,
0077         x = randperm(N);
0078         group = group(x);
0079         for j = 1:r,
0080             mu(j) = CircularMean(angles(group==j));
0081         end
0082         fr_rand(i) = ComputeStatistics(angles,group,mu);
0083     end
0084     % Sort statistics
0085     fr_rand = sort(fr_rand);
0086     % Where does the statistics for the original data fall?
0087     m = find(fr_rand>=fr);
0088     % p level
0089     if isempty(m),
0090         p = 0;
0091     else
0092         m0 = sum(m==m(1));
0093         m = m(1);
0094         if m0 == 1,
0095             p = (nRandomizations-m+1)/nRandomizations;
0096         else
0097             p = (nRandomizations-m)/nRandomizations+m0/2;
0098         end
0099     end
0100 end
0101 
0102 h = p < alpha;
0103 
0104 message = ['Concentration test: median kappa = ' num2str(kappa) ' ('];
0105 for i = 1:r,
0106     message = [message num2str(k(i)) ' '];
0107 end
0108 message(end) = ')';
0109 disp(message);
0110 
0111 if h,
0112     message = '+++ Two or more concentration parameters are significantly different';
0113 else
0114     message = '--- Concentration parameters are not significantly different';
0115 end
0116 message = ['Concentration test: ' message ' (p='  num2str(p) ', fr=' num2str(fr) ', N=' int2str(N) ')'];
0117 disp(message);
0118 
0119 
0120 function fr = ComputeStatistics(angles,group,mu)
0121 
0122 global n N r;
0123 
0124 for i = 1:r,
0125     a = abs(sin(angles(group==i)-mu(i)));
0126     d(i) = sum(a)/n(i);
0127     s(i) = sum((a-d(i)).^2);
0128 end
0129 
0130 D = sum(n.*d)/N;
0131 
0132 fr = (N-r)*sum(n.*(d-D).^2)/((r-1)*sum(s));