CircularConfidenceIntervals

CircularConfidenceIntervals - Compute circular mean and confidence intervals.

 Compute circular mean and confidence intervals of circular data (see Fisher,
 Analysis of circular data, p. 88). If fewer than 25 angle measurements are
 provided, bootstrap is used.

  USAGE

    [mean,boundaries] = CircularConfidenceIntervals(angles,alpha,nBootstrap)

    angles         angles in radians
    alpha          optional confidence level (default = 0.05)
    nBootstrap     optional number of bootstraps (for < 25 angle values)
                   (default = 2000)

  SEE

    See also CircularMean, CircularVariance, Concentration, ConcentrationTest.

0001 function [m,boundaries] = CircularConfidenceIntervals(angles,alpha,nBootstrap)
0002 
0003 %CircularConfidenceIntervals - Compute circular mean and confidence intervals.
0004 %
0005 % Compute circular mean and confidence intervals of circular data (see Fisher,
0006 % Analysis of circular data, p. 88). If fewer than 25 angle measurements are
0007 % provided, bootstrap is used.
0008 %
0009 %  USAGE
0010 %
0011 %    [mean,boundaries] = CircularConfidenceIntervals(angles,alpha,nBootstrap)
0012 %
0013 %    angles         angles in radians
0014 %    alpha          optional confidence level (default = 0.05)
0015 %    nBootstrap     optional number of bootstraps (for < 25 angle values)
0016 %                   (default = 2000)
0017 %
0018 %  SEE
0019 %
0020 %    See also CircularMean, CircularVariance, Concentration, ConcentrationTest.
0021 
0022 % Copyright (C) 2004-2014 by Michaël Zugaro
0023 %
0024 % This program is free software; you can redistribute it and/or modify
0025 % it under the terms of the GNU General Public License as published by
0026 % the Free Software Foundation; either version 3 of the License, or
0027 % (at your option) any later version.
0028 
0029 if nargin < 1,
0030     error('Incorrect number of parameters (type ''help <a href="matlab:help CircularConfidenceIntervals">CircularConfidenceIntervals</a>'' for details).');
0031 end
0032 isradians(angles);
0033 if isdscalar(angles),
0034     m = angles;
0035     boundaries = [m m];
0036     return
0037 end
0038 
0039 m = [];
0040 boundaries = [];
0041 n = length(angles);
0042 if n == 0, return; end
0043 
0044 if nargin < 2
0045     alpha = 0.05;
0046 end
0047 if nargin < 3
0048     if n < 25
0049         nBootstrap = 2000;
0050     else
0051         nBootstrap = 0;
0052     end
0053 end
0054 
0055 if nBootstrap == 0
0056     [m,r1] = CircularMean(angles);
0057     [~,r2] = CircularMean(wrap(2*angles));
0058     delta = (1-r2)./(2*r1.^2);
0059     sigma = sqrt(delta/n);
0060     sinarg = norminv(1-alpha/2)*sigma;
0061     if sinarg < 1
0062         err = asin(sinarg);
0063     else
0064         err = pi;
0065     end
0066     boundaries = [m-err; m+err];
0067 else
0068     m = CircularMean(angles);
0069     b = bootstrp(nBootstrap,'CircularMean',angles);
0070 
0071     % unwrap data around mean value
0072     unwrapped = mod(b-repmat(m,nBootstrap,1)+pi,2*pi)+repmat(m,nBootstrap,1)-pi;
0073 
0074     boundaries(1,:) = prctile(unwrapped,100*(alpha/2));
0075     boundaries(2,:) = prctile(unwrapped,100-100*(alpha/2));
0076 end