Home > atmlab > circular > circ_kappa.m

circ_kappa

PURPOSE ^

SYNOPSIS ^

function kappa = circ_kappa(alpha,w)

DESCRIPTION ^

 kappa = circ_kappa(alpha,[w])
   Computes an approximation to the ML estimate of the concentration 
   parameter kappa of the von Mises distribution.

   Input:
     alpha   angles in radians OR alpha is length resultant
     [w      number of incidences in case of binned angle data]

   Output:
     kappa   estimated value of kappa

   References:
     Statistical analysis of circular data, Fisher, equation p. 88

 Circular Statistics Toolbox for Matlab

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

DOWNLOAD ^

circ_kappa.m

SOURCE CODE ^

0001 function kappa = circ_kappa(alpha,w)
0002 %
0003 % kappa = circ_kappa(alpha,[w])
0004 %   Computes an approximation to the ML estimate of the concentration
0005 %   parameter kappa of the von Mises distribution.
0006 %
0007 %   Input:
0008 %     alpha   angles in radians OR alpha is length resultant
0009 %     [w      number of incidences in case of binned angle data]
0010 %
0011 %   Output:
0012 %     kappa   estimated value of kappa
0013 %
0014 %   References:
0015 %     Statistical analysis of circular data, Fisher, equation p. 88
0016 %
0017 % Circular Statistics Toolbox for Matlab
0018 
0019 % By Philipp Berens, 2009
0020 % berens@tuebingen.mpg.de - www.kyb.mpg.de/~berens/circStat.html
0021 
0022 
0023 alpha = alpha(:);
0024 
0025 if nargin<2
0026   % if no specific weighting has been specified
0027   % assume no binning has taken place
0028     w = ones(size(alpha));
0029 else
0030   if size(w,2) > size(w,1)
0031     w = w';
0032   end 
0033 end
0034 
0035 N = length(alpha);
0036 
0037 if N>1
0038   R = circ_r(alpha,w);
0039 else
0040   R = alpha;
0041 end
0042 
0043 if R < 0.53
0044   kappa = 2*R + R^3 + 5*R^5/6;
0045 elseif R>=0.53 && R<0.85
0046   kappa = -.4 + 1.39*R + 0.43/(1-R);
0047 else
0048   kappa = 1/(R^3 - 4*R^2 + 3*R);
0049 end
0050 
0051 if N<15 && N>1
0052   if kappa < 2
0053     kappa = max(kappa-2*(N*kappa)^-1,0);    
0054   else
0055     kappa = (N-1)^3*kappa/(N^3+N);
0056   end
0057 end

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