COO-207 Histogram with Steve Watts' Shannon entropy based bin counting algorithm

[prkrtg]

Paper: https://arxiv.org/pdf/2210.02848

[weatherhead99]

couple of minor things to look at and reply to. implementation looks of good quality
and nice test coverage. No major objections from me.

[weatherhead99]

nice, was going to suggest something like this from the conversation we had yesterday about the narrowing/ / widening float copy

[weatherhead99]

forgive me if I'm missing something here but why do we need counts to be floating point?

In fact, don't we specifically want it NOT to be? How do you have half a count of something

[prkrtg]

Yea good catch, it's now counts = np.asarray(counts)

[weatherhead99]

here I think the data may well be required to be floating point. We could for example use

def differential_entropy_knn(data: npt.NDArray[np.floating]) : ...

which will allow mypy to catch a more specific case

[weatherhead99]

I know you won't like this, but maybe we should consider not implementing KNN ourselves and instead leveraging an entropy based KNN from scikit-learn https://scikit-learn.org/stable/modules/generated/sklearn.feature_selection.mutual_info_regression.html

Or maybe that's way overkill. I find it very likely we'll eventually end up depending on scikit-image, but perhaps not scikit-learn. That said it's not all that complicated and maybe that is overkill

[prkrtg]

I scoped out sklearn.mutual_info_regression , the problem is it estimates information I(X;Y) between two variables, not the marginal differential entropy h(X) of a single distribution from Eq. 11. Basically we can't just use this function to do I(X;Y) = h(X) + h(Y) - h(X,Y) back into just h(X)

So scipy.stats.differential_entropy does implement it, it passes the basic tests

import time
import numpy as np
import scipy.stats as st
import warnings

rng = np.random.default_rng(0)

# performance for 8k x 8k scale
n = 8192 * 8192
data = rng.normal(1000, 5, n)
t0 = time.perf_counter()
h = st.differential_entropy(data, method='auto')
t1 = time.perf_counter()
print(f'scipy differential_entropy on {n:,} points: {t1-t0:.2f}s, h={h/np.log(2):.4f} bits')

# behavior on quantized/tied data
data_q = np.round(rng.normal(1000, 5, 100_000))
print('n exact ties in quantized sample:', (np.diff(np.sort(data_q)) == 0).sum(), 'out of', data_q.size)
with warnings.catch_warnings(record=True) as w:
    warnings.simplefilter('always')
    h_q = st.differential_entropy(data_q, method='auto')
    for warning in w:
        print('WARNING:', warning.category.__name__, warning.message)
print('scipy result on quantized (tied) data:', h_q / np.log(2), 'bits')

This spits out -inf and also RuntimeWarning: divide by zero, which makes sense since as a number gets smaller and smaller toward zero, its log heads toward negative infinity. So I thought since real detector images are quantized, there are tons of exact-zero gaps which is why in the home-made knn function I have the dither option.

But let me know your thoughts

[weatherhead99]

consider a dataclass rather than a dict here, it's trivially convertible if needed and semantically preferable (at least IMO)

[weatherhead99]

if dict has to be used, please change type annotation to ->dict[str, Any]

[prkrtg]

Long term I think a data class would be good, but looking at the other areas of eregion like image_stats.py which currently all return plain dicts, I'll keep it as a dict
For now changed type annotation to ->dict[str, Any]

[weatherhead99]

might be worth using a pytest harness to set up the RNG, saves a few lines.
Might not be, also

[prkrtg]

Idk if it's worth it (at least not yet), since I would only really replace np.random.default_rng(N), I still have to dp sample = rng.normal or .uniform, etc with different distributions and sizes

Maybe if this repeated pattern is in more tests