Add: Angular Randomisation Test

[elhananby]

An implementation of the Angular Randomisation test (ART), which was originally introduced by Ali & Abushilah (2022) and later validated by Ruxton, Malkemper & Landler (2023). Their validation found that ART generally outperforms traditional two-sample tests like Watson's U² and Watson-Wheeler tests, although not for all use cases.

I've added a cell to the hypothesis notebook, although I have not implemented a proper test yet.

References:

Ali, A. J., & Abushilah, S. F. (2022). Distribution-free two-sample homogeneity test for circular data based on geodesic distance. Int. J. Nonlinear Anal. Appl., 13, 2703-2711.

Ruxton, G. D., Malkemper, E. P., & Landler, L. (2023). Evaluating the power of a recent method for comparing two circular distributions: an alternative to the Watson U² test. Scientific Reports, 13, 10007.

[huangziwei]

thanks! will have a look soon.

[huangziwei]

@elhananby Thanks a lot! This PR is a great addition. I haven't gotten around to reading the papers yet, but it looks good. I left a few comments—if you'd like, you can make the changes, or I can do it next week with a PR to your fork, and then we can merge it.

P.S. It would be great if you could pip install watermark and rerun the notebook. That will mark the versions of the packages you used (it currently raised an error because you don’t have it installed). Also, please place the new cell above the watermark. :)

[huangziwei]

Love this. I am actually thinking about adding more distance metrics to the descriptive.py. Geodesic distance is nice to have, along with other metrics in dist.circular(). Here's a draft for the new circ_dist and circ_pairdist:

def circ_dist(x: np.ndarray, y: Optional[np.ndarray] = None, metric: str = "center", return_sum: bool = False) -> np.ndarray:
    r"""
    Compute the element-wise circular distance between two arrays of angles.

    Parameters
    ----------
    x : array-like
        First sample of circular data (radians).
    y : array-like, optional
        Second sample of circular data (radians). If None, computes element-wise 
        distances within `x` itself.
    metric : str, optional
        Distance metric to use, options:
        - "center" (default): Standard circular difference wrapped to [-π, π].
        - "geodesic": π - |π - |x - y||.
        - "angularseparation": 1 - cos(x - y).
        - "chord": sqrt(2 * (1 - cos(x - y))).
    return_sum : bool, optional
        If True, returns the sum of all computed distances (like R's `dist.circular()`).

    Returns
    -------
    array
        Element-wise distance values based on the chosen metric.
    """
    x = np.asarray(x)

    if y is None:
        y = x

    y = np.asarray(y)

    # Ensure broadcasting works without explicit shape checks
    try:
        np.broadcast_shapes(x.shape, y.shape)
    except ValueError:
        raise ValueError(f"Shapes {x.shape} and {y.shape} are incompatible for broadcasting.")

    if metric == "center":
        distances = np.angle(np.exp(1j * x) / np.exp(1j * y))

    elif metric == "geodesic":
        distances = np.pi - np.abs(np.pi - np.abs(x - y))

    elif metric == "angularseparation":
        distances = 1 - np.cos(x - y)

    elif metric == "chord":
        distances = np.sqrt(2 * (1 - np.cos(x - y)))

    else:
        raise ValueError(f"Unknown metric: {metric}")
    
    return np.sum(distances) if return_sum else distances


def circ_pairdist(x: np.ndarray, y: Optional[np.ndarray] = None, metric: str = "center", return_sum: bool = False) -> np.ndarray:
    r"""
    Compute the pairwise circular distance between all elements in `x` and `y`.

    Parameters
    ----------
    x : array-like
        First sample of circular data (radians).
    y : array-like, optional
        Second sample of circular data (radians). If None, computes pairwise 
        distances within `x` itself.
    metric : str, optional
        Distance metric to use (same options as `circ_dist`).
    return_sum : bool, optional
        If True, returns the sum of all computed distances (like R's `dist.circular()`).

    Returns
    -------
    ndarray
        Pairwise distance matrix where entry (i, j) is the circular distance 
        between x[i] and y[j] based on the chosen metric.
    """
    x = np.asarray(x)

    # If y is not provided, compute pairwise distances within x
    if y is None:
        y = x

    y = np.asarray(y)

    # Reshape to allow broadcasting for pairwise computation
    x_reshaped = x[:, None]  # Shape (n, 1)
    y_reshaped = y[None, :]  # Shape (1, m)

    return circ_dist(x_reshaped, y_reshaped, metric=metric, return_sum=return_sum)

and then you can remove this inner function and replace the line total_distance = compute_geodesic_distance(S1, S2) below with:

# from pycircstat2.descriptive import circ_pairdist
total_distance = circ_pairdist(S1, S2, metric="geodesic", return_sum=True)

do you want to add it for me or I should make a pr to your fork first?

[elhananby]

Hi @huangziwei, I've made the changes and pushed them to my fork.

[huangziwei]

I haven't gotten around to reading the papers yet, but a simple test for input/output to make sure the code runs would be enough.

[elhananby]

done in latest push.

[huangziwei]

Looks good already. I will just merge it now. Thanks! :)