Merge pull request #3 from dschmitz89/1.2

1.2

README.md

@@ -19,7 +19,7 @@ Refer to the [online documentation](https://spherical-stats.readthedocs.io/en/la
     * Orientation tensor
 * Parametric distributions with scipy.stats like API:
     * Modeling axial data: Angular central gaussian distribution (ACG)
-    * Modeling vector data: Elliptically symmetrical angular gausian distribution (ESAG)
+    * Modeling vector data: Elliptically symmetrical angular gausian distribution (ESAG), Von Mises-Fisher distribution (VMF)
 
 Example usage of the distributions:
 

setup.py

@@ -7,7 +7,7 @@
 
 setup(
     name='spherical_stats',   
-    version='1.1', 
+    version='1.2', 
     description='Spherical statistics in Python',
     author='Daniel Schmitz',
     author_email='[email protected]',

spherical_stats/_acg.py

@@ -66,12 +66,21 @@ def fit(vectors, tol):
     return l
 
 class ACG(object):
-    '''
+    r'''
     Angular Central Gaussian distribution
 
     Args:
         cov_matrix (optional, ndarray (3, 3) ): Covariance matrix of the ACG distribution
 
+    The angular central gaussian distribution is an elliptically symmetrical distribution
+    for axial data. Its PDF is defined as 
+
+    .. math::
+
+        p_{ACG}(\mathbf{x}|\mathbf{\Lambda}) = \frac{1}{4\pi\sqrt{|\Lambda|}}(\mathbf{x}\Lambda^{-1}\mathbf{x})^{-\frac{3}{2}}
+
+    with covariance matrix :math:`\Lambda` and its determinant :math:`|\Lambda|` .
+
     Notes
     -------
     Reference: Tyler, Statistical  analysis  for  the  angular  central Gaussian  distribution  on  the  sphere

spherical_stats/_esag.py

@@ -222,22 +222,27 @@ def _fit(vectors, print_summary = False):
     return optimized_params
 
 class ESAG(object):
-    '''
+    r'''
     Elliptically symmetrical angular Central Gaussian distribution
 
     Args:
         params (optional, ndarray (5, ) ): Parameters of the distribution
         
-    ``params`` are the following: :math:`(\mu_0, \mu_1, \mu_2, \gamma_1, \gamma_2)`.\n    
-    The principal orientation vectors is given by the normalized vector :math:`\mu=(\mu_0, \mu_1, \mu_2)/||(\mu_0, \mu_1, \mu_2)||`.\n
+    ``params`` are the following: :math:`(\mu_0, \mu_1, \mu_2, \gamma_1, \gamma_2)`.
+
+    The principal orientation vectors is given by the normalized vector :math:`\boldsymbol{\mu}=[\mu_0, \mu_1, \mu_2]^T/||[\mu_0, \mu_1, \mu_2]^T||` 
     and the shape of the distribution is controlled by the parameters :math:`\gamma_1` and :math:`\gamma_2`.
 
     Notes
     -------
-    Reference: Paine et al. An elliptically symmetric angular Gaussian distribution, 
-    Statistics and Computing volume 28, 689–697 (2018)\n
-    Unlike the original Matlab implementation the distribution is fitted using the L-BFGS-B algorithm
+    The formula of the ESAG PDF is quite complicated, developers are referred to the reference below. 
+
+    Note that unlike the original Matlab implementation the distribution is fitted using the L-BFGS-B algorithm
     based on a finite difference approximation of the gradient. So far, this has proven to work succesfully.
+
+    Reference: Paine et al. An elliptically symmetric angular Gaussian distribution, 
+    Statistics and Computing volume 28, 689–697 (2018)
+
     '''
 
     def __init__(self, params = None):

spherical_stats/_vmf.py

@@ -67,19 +67,32 @@ def obj(s):
     return mu, kappa
 
 class VMF:
-    '''
-    Van Mises-Fisher distribution
+    r"""
+    Von Mises-Fisher distribution
 
     Args:
-        mu (optional, ndarray (3, ) ): Mean orientation
+        mu (optional, ndarray (3, ) ): Mean orientation 
         kappa (optional, float): positive concentration parameter
 
-    Notes
-    -------
-    References:\n
-    Mardia, Jupp. Directional Statistics, 1999. \n
-    Numerically stable sampling of the von Mises Fisher distribution on  S2. Wenzel, 2012
-    '''
+    The VMF distribution is an isotropic distribution for 
+    directional data. Its PDF is typically defined as 
+
+    .. math::
+
+        p_{vMF}(\mathbf{x}| \boldsymbol{\mu}, \kappa) = \frac{\kappa}{4\pi\cdot\text{sinh}(\kappa)}\exp(\kappa \boldsymbol{\mu}^T\mathbf{x})
+
+    Here, the numerically stable variant from (Wenzel, 2012) is used:
+
+    .. math::
+
+        p_{vMF}(\mathbf{x}| \boldsymbol{\mu}, \kappa) = \frac{\kappa}{2\pi(1-\exp(-2\kappa))}\exp(\kappa( \boldsymbol{\mu}^T \mathbf{x}-1))
+    
+    References:
+
+    Mardia, Jupp. Directional Statistics, 1999. 
+
+    Wenzel. Numerically stable sampling of the von Mises Fisher distribution on  S2. 2012
+    """
     def __init__(self, mu = None, kappa = None):
 
         self.mu = mu
@@ -146,4 +159,40 @@ def fit(self, data):
         data : ndarray (n, 3)
             Vector data the distribution is fitted to
         '''
-        self.mu, self.kappa = _fit(data)    
+        self.mu, self.kappa = _fit(data)
+
+    def angle_within_probability_mass(self, alpha, deg = False):
+        r"""
+        Calculates the angle which contains probability mass alpha of the VMF density around the mean angle
+
+        Reference: Fayat, 2021. Conversion of the von Mises-Fisher concentration parameter to an equivalent angle.
+
+        https://github.com/rfayat/SphereProba/blob/main/ressources/vmf_integration.pdf
+
+        Arguments
+        ----------
+        alpha : float
+            Probability mass. Must be :math:`0<\alpha<1` 
+        deg : optional, default False
+            If True, converts the result into degrees
+
+        Returns
+        ----------
+        angle : float
+            Resulting angle
+        
+        """
+        if self.kappa is not None:
+
+            nominator = np.log(1-alpha + alpha * np.exp(-2 * self.kappa))  
+
+            angle = np.arccos(1+nominator/self.kappa)
+
+            if deg == True:
+                angle=np.rad2deg(angle)
+
+            return angle
+
+        else:
+
+            raise ValueError("Concentration parameter kappa unknown.")