Probability density estimation in circular variables

[gerompampastrumf]

Description
The current continuous estimators (e.g., KSG, kernel-based methods) assume that variables lie in a Euclidean space (ℝⁿ). However, in many practical applications, variables are periodic rather than linear. A common example is the phase of oscillatory signals, which wraps at 2π (or equivalently at −π and π). In such cases, the underlying geometry is circular (S¹)

Suggestion
It would be nice if there could be a distance_space argument in the estimators, which would allow for Euclidean or circular spaces.

References where information theory metrics are used on phases
Phase Transfer Entropy: https://doi.org/10.1016/j.neuroimage.2013.08.056
Phase Mutual Information: https://doi.org/10.1371/journal.pone.0044633

[cbueth]

Thanks for this suggestion and references. The extension for periodic boundary conditions sounds intuitive. The question would be, where does this fit into the package structure, so that most of the estimators can profit off of it. Feel free to open a PR suggesting an implementation. We can also have a chat to find out where to best place it.

[kjohnsen]

More useful than being able to set a single distance_space argument would be to allow for different spaces/kernels per variable. In my case, for example, all my variables are (not circular) continuous except for one of my conditioning variables, which is discrete. I'd like to be able to combine the typical Gaussian kernel for the continuous variables with the Aitchison-Aitken kernel or the Dirac delta kernel for the discrete variable

[cbueth]

Do you know how these spaces are usually combined (with any estimation technique), is there an application which sets this up or/and paper defining the mixture of spaces? Personally, I have not worked with other spaces.
Just assuming, for some input spaces this might come down simple data preparation, i.e. to a transformation like a modulo, before passing to the existing estimators. This would be the simple case, but I'm unsure which are the pitfalls to watch out to. I'm happy to hear your input so we might nail the requirements down here for a solution serving multiple usecases.

[kjohnsen]

I am by no means an expert on this topic so I can't point to a paper---I've stumbled my way to this point with AI conversations. The reason I want to use kernel estimation is because I have low sample sizes. From what I understand, I could treat my discrete conditioning variable as continuous, but that would be kind of hacky (like, why is random seed 1 closer to 2 than to 3?). Or I could compute MI separately for each condition and average together, but that would be very noisy when I have small per-condition group sizes (some as small as 2-4 samples). Aitchison-Aitken allows me to smooth together info from different groups a little bit but still keep them mostly separate.

From what I understand, combining continuous kernels with discrete kernels would look like this
$$p(x_i, y_i, z_i) = \frac{1}{n} \sum_{j=1}^{n} \left[ K_{cont}\left(dist_{cont}(i, j)\right) \times \prod_{z_{disc}} K_{disc}(z_{disc,i}, z_{disc,j}) \right]$$

This could get really complicated if you want maximal flexibility since in some cases you might not even want to group all your continuous variables: when you have a mix of circular and linear, when you want to use different kernels and/or bandwidths for different variables/groups of variables, etc. Not sure how to implement this with a clean, user-friendly interface. Probably want to keep the easy-to-use current behavior of specifying a single kernel/bandwidth for all, but optionally accept kernel groups, e.g., [...('kernel_type', bandwidth, [indices of vars in this group])]

For example

# the one-kernel-fits-all default:
im.cmi(x, y, cond=z, approach='kernel', kernel='box', bandwidth=1.5)

# the separate kernel groups idea
im.cmi(x_space1, x_space2, x_circular, y, cond=(z1, z2), approach="kernel", kernel_groups=[
    ('gaussian', 0.7, [0,1]), 
    ('von mises', 'scott', [2]),
    ('gaussian', 0.4, [3]),
    # for discrete Z, though X's could be discrete or Z's could be continuous too
    ('AA', 0.1, [4]),
    ('dirac', None, [5]),
])

Oh, and I threw in the idea of bandwidth accepting strings like 'scott' to estimate the optimal bandwidth

[cbueth]

Thanks for sharing these thoughts. I see there might be quite some applications and having one usable interface might be quite neet. Right now I do not have any better one than the one you suggest. It could be a strength if one separates the data and estimation/coincidence counting of marginal spaces before calculating joint measures. So, an alternative would be let the user create data + estimation pairs before passing it to the measure. But honestly, I would also keep the current function interfaces, as they largely simplify most use cases and have been well received.
One recent package that implements mixed measures is this one: https://syntropy.readthedocs.io/en/latest/api/syntropy.mixed.html. I think it really come down on how it is best implemented internally in infomeasure, keeping efficiency high and code maintenance low.

The coming month I am not able to work on this doe to other tasks, but am happy to review PRs or to give feedback and answer. @kjohnsen, it would be nice if you could open a separate issue as freature request for this, so we can keep this issue by @gerompampastrumf separated, while somewhat connected, to it.

[gerompampastrumf]

Returning back to the problem of estimation in circular variables, I have explored the code and I believe that there are two distinct cases to be addressed: a simple one and a more general one.

1. In the simple case, all variables are phases bounded within the periodic interval ([0, L]). Then a minimal change to the scipy.spatial.KDTree function used for the nearest neighbor search is enough to address this case. More specifically, we only need to pass boxsize=L as an argument in all calls of KDTree. For example,
   tree_joint = KDTree(joint_space_data)
   should become
   tree_joint = KDTree(joint_space_data, boxsize=self.boxsize).

2. In a more general case, it would make sense to have the freedom to include some continuous variables with periodic boundary conditions and others with non-periodic conditions. While KDTree allows for different box sizes, it does not allow for a mixture of periodic and non-periodic variables. Nonetheless, we can utilize the fact that if we pass a sufficiently large value in the boxsize argument (e.g., 2*max(x) where x is a variable in [0, max(x)]), then these variables are effectively treated as non-periodic.

I have modified the kraskov_stoegbauer_grassberger estimator of TE according to point 1 locally. However, I have not yet made a pull request because I have not tested it thoroughly, and it is currently only applied to a single estimator. Also, if we decide to properly add this feature, it would make sense to proceed directly with point 2.

[cbueth]

That sounds good. I think this addition only makes sense for continuous data, so it could be added to these estimators directly (i.e. kNN, KL/KSG). In this case, the addition to KDTree and adding a new keyword argument is the way to go, and @kjohnsen feature request can be handled separately. When you have the PR done, I'll look into it.