Home > atmlab > circular > circ_raotest.m

circ_raotest

PURPOSE ^

[p U UC] = circ_raotest(alpha)

SYNOPSIS ^

function [p U UC] = circ_raotest(alpha)

DESCRIPTION ^

 [p U UC] = circ_raotest(alpha)
   Calculates Rao's spacing test by comparing distances between points on
   a circle to those expected from a uniform distribution.

   H0: Data is distributed uniformly around the circle. 
   H1: Data is not uniformly distributed around the circle.

   Alternative to the Rayleigh test and the Omnibus test. Less powerful
   than the Rayleigh test when the distribution is unimodal on a global
   scale but uniform locally.

   Due to the complexity of the distributioin of the test statistic, we
   resort to the tables published by 
       Russell, Gerald S. and Levitin, Daniel J.(1995)
       'An expanded table of probability values for rao's spacing test'
       Communications in Statistics - Simulation and Computation
   Therefore the reported p-value is the smallest alpha level at which the
   test would still be significant. If the test is not significant at the
   alpha=0.1 level, we return the critical value for alpha = 0.05 and p =
   0.5.

   Input:
     alpha     sample of angles

   Output:
     p         smallest p-value at which test would be significant
     U         computed value of the test-statistic u
     UC        critical value of the test statistic at sig-level


   References:
     Batschelet, 1981, Sec 4.6

 Circular Statistics Toolbox for Matlab

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SUBFUNCTIONS ^

DOWNLOAD ^

circ_raotest.m

SOURCE CODE ^

0001 function [p U UC] = circ_raotest(alpha)
0002 
0003 % [p U UC] = circ_raotest(alpha)
0004 %   Calculates Rao's spacing test by comparing distances between points on
0005 %   a circle to those expected from a uniform distribution.
0006 %
0007 %   H0: Data is distributed uniformly around the circle.
0008 %   H1: Data is not uniformly distributed around the circle.
0009 %
0010 %   Alternative to the Rayleigh test and the Omnibus test. Less powerful
0011 %   than the Rayleigh test when the distribution is unimodal on a global
0012 %   scale but uniform locally.
0013 %
0014 %   Due to the complexity of the distributioin of the test statistic, we
0015 %   resort to the tables published by
0016 %       Russell, Gerald S. and Levitin, Daniel J.(1995)
0017 %       'An expanded table of probability values for rao's spacing test'
0018 %       Communications in Statistics - Simulation and Computation
0019 %   Therefore the reported p-value is the smallest alpha level at which the
0020 %   test would still be significant. If the test is not significant at the
0021 %   alpha=0.1 level, we return the critical value for alpha = 0.05 and p =
0022 %   0.5.
0023 %
0024 %   Input:
0025 %     alpha     sample of angles
0026 %
0027 %   Output:
0028 %     p         smallest p-value at which test would be significant
0029 %     U         computed value of the test-statistic u
0030 %     UC        critical value of the test statistic at sig-level
0031 %
0032 %
0033 %   References:
0034 %     Batschelet, 1981, Sec 4.6
0035 %
0036 % Circular Statistics Toolbox for Matlab
0037 
0038 % By Philipp Berens, 2009
0039 % berens@tuebingen.mpg.de
0040 
0041 alpha = alpha(:);
0042 
0043 % for the purpose of the test, convert to angles
0044 alpha = circ_rad2ang(alpha);
0045 n = length(alpha);
0046 alpha = sort(alpha);
0047 
0048 % compute test statistic
0049 U = 0;
0050 lambda = 360/n;
0051 for j = 1:n-1
0052     ti = alpha(j+1) - alpha(j);
0053     U = U + abs(ti - lambda);
0054 end
0055 
0056 tn = (360 - alpha(n) + alpha(1));
0057 U = U + abs(tn-lambda);
0058 
0059 U = (1/2)*U;
0060 
0061 % get critical value from table
0062 [p UC] = getVal(n,U);
0063 
0064 
0065 
0066 function [p UC] = getVal(N, U)
0067 
0068 % Table II from Russel and Levitin, 1995
0069 
0070 alpha = [0.001, .01, .05, .10];
0071 table = [ 4   247.32, 221.14, 186.45, 168.02;
0072           5   245.19, 211.93, 183.44, 168.66;
0073           6   236.81, 206.79, 180.65, 166.30;
0074           7   229.46, 202.55, 177.83, 165.05;
0075           8   224.41, 198.46, 175.68, 163.56;
0076           9   219.52, 195.27, 173.68, 162.36;
0077           10  215.44, 192.37, 171.98, 161.23;
0078           11  211.87, 189.88, 170.45, 160.24;
0079           12  208.69, 187.66, 169.09, 159.33;
0080           13  205.87, 185.68, 167.87, 158.50;
0081           14  203.33, 183.90, 166.76, 157.75;
0082           15  201.04, 182.28, 165.75, 157.06;
0083           16  198.96, 180.81, 164.83, 156.43;
0084           17  197.05, 179.46, 163.98, 155.84;
0085           18  195.29, 178.22, 163.20, 155.29;
0086           19  193.67, 177.08, 162.47, 154.78;
0087           20  192.17, 176.01, 161.79, 154.31;
0088           21  190.78, 175.02, 161.16, 153.86;
0089           22  189.47, 174.10, 160.56, 153.44;
0090           23  188.25, 173.23, 160.01, 153.05;
0091           24  187.11, 172.41, 159.48, 152.68;
0092           25  186.03, 171.64, 158.99, 152.32;
0093           26  185.01, 170.92, 158.52, 151.99;
0094           27  184.05, 170.23, 158.07, 151.67;
0095           28  183.14, 169.58, 157.65, 151.37;
0096           29  182.28, 168.96, 157.25, 151.08;
0097           30  181.45, 168.38, 156.87, 150.80;
0098           35  177.88, 165.81, 155.19, 149.59;
0099           40  174.99, 163.73, 153.82, 148.60;
0100           45  172.58, 162.00, 152.68, 147.76;
0101           50  170.54, 160.53, 151.70, 147.05;
0102           75  163.60, 155.49, 148.34, 144.56;
0103           100 159.45, 152.46, 146.29, 143.03;
0104           150 154.51, 148.84, 143.83, 141.18;
0105           200 151.56, 146.67, 142.35, 140.06;
0106           300 148.06, 144.09, 140.57, 138.71;
0107           400 145.96, 142.54, 139.50, 137.89;
0108           500 144.54, 141.48, 138.77, 137.33;
0109           600 143.48, 140.70, 138.23, 136.91;
0110           700 142.66, 140.09, 137.80, 136.59;
0111           800 142.00, 139.60, 137.46, 136.33;
0112           900 141.45, 139.19, 137.18, 136.11;
0113           1000  140.99, 138.84, 136.94, 135.92  ];
0114         
0115 ridx = find(table(:,1)>=N,1);    
0116 cidx = find(table(ridx,2:end)<U,1);
0117 
0118 if ~isempty(cidx)
0119   UC = table(ridx,cidx+1);
0120   p = alpha(cidx);
0121 else
0122   UC = table(ridx,end-1);
0123   p = .5;
0124 end
0125 
0126 
0127 
0128 
0129 
0130

Generated on Mon 15-Sep-2014 13:31:28 by m2html © 2005