Use circstd**2 for variance to ensure consistency across SciPy versions

[ale94mleon]

Problem

In SciPy >=1.8, the circvar function was changed to return a normalized variance (bounded [0,1]), whereas older versions returned an unnormalized variance. This breaks the _sum_var() calculation in the repository, because bond (linear) variances are still absolute, and summing them with normalized angular variances no longer makes sense.

How to reproduce

Installing different SciPy versions

micromamba create -n scipy-less-1.8 'scipy<1.8' # installed scipy=1.7.3
micromamba create -n scipy-newest scipy # installed scipy=1.16.1

Testing script

The following script shows the difference in behavior between SciPy <1.8
and SciPy >=1.8:

import warnings
with warnings.catch_warnings():
    warnings.simplefilter("ignore", RuntimeWarning)
    import numpy as np
    from scipy.stats import circmean, circstd, circvar
    from scipy import __version__ as scipy_version

# Sample angles in degrees
angles = np.array([10, 20, 30, 350])

# -----------------------------
# Circular statistics
# -----------------------------
circmean_v = circmean(angles, high=360, low=0)
circstd_v = circstd(angles, high=360, low=0)
circvar_v = circvar(angles, high=360, low=0)
circstd2_v = circstd_v**2

abs_diff = abs(circstd2_v - circvar_v)
rel_diff = abs_diff / (abs(circvar_v) + 1e-16)  # avoid divide by zero
similar = np.isclose(circstd2_v, circvar_v, rtol=1e-8, atol=1e-10)

print(f"\n\nSciPy version: {scipy_version}")
print("=== Circular Statistics ===")
print(f"circmean: {circmean_v}")
print(f"circstd: {circstd_v}")
print(f"circvar: {circvar_v}")
print(f"circstd^2: {circstd2_v}")
print(f"Absolute diff (circstd^2 - circvar): {abs_diff}")
print(f"Relative diff: {rel_diff}")
print(f"Are circstd^2 and circvar similar? {similar}\n")

# -----------------------------
# Linear statistics
# -----------------------------
mean_v = angles.mean()
std_v = angles.std()
var_v = angles.var()
std2_v = std_v**2

abs_diff_lin = abs(std2_v - var_v)
rel_diff_lin = abs_diff_lin / (abs(var_v) + 1e-16)
similar_lin = np.isclose(std2_v, var_v, rtol=1e-8, atol=1e-10)

print("=== Linear Statistics ===")
print(f"mean: {mean_v}")
print(f"std: {std_v}")
print(f"var: {var_v}")
print(f"std^2: {std2_v}")
print(f"Absolute diff (std^2 - var): {abs_diff_lin}")
print(f"Relative diff: {rel_diff_lin}")
print(f"Are std^2 and var similar? {similar_lin}")

Running script

micromamba run -n 'scipy-less-1.8' python script.py
micromamba run -n 'scipy-newest' python script.py

Results

• In SciPy=1.7.3, circvar ≈ circstd**2.
• In SciPy =1.16.1, circvar is normalized (~0.03 in the example), so combining it with bond variances gives incorrect results.

Solution

This PR updates the VectorData.analyze() method to calculate variance as circstd**2 for periodic data:

self.var = self.stdev**2

And similarly for linear data for consistency.

This ensures _sum_var() produces the expected combined variance consistently across all supported SciPy versions.

Related issue

Fixes #41

⚠️ Observation

The terms in _sum_var() (for bonds and angles) should ideally be in the same units. Currently, it seems that angles are in degrees rather than radians. Not sure if this is intentional. Therefore an extra conversion step might be needed.

[IAlibay]

Thanks for the PR @ale94mleon - looks like there's some issues with the package build that need to be handled before we merge this. I'll work on them in the next couple of days.

[IAlibay]

Thanks, this seems good to me, apologies for the delay.

I'll merge this and deal with the package build issues afterwards, it'll make life easier.