Probability: tfp.stats.percentile sort order differs from docstring (and differs in behavior from NumPy)

[matthewfeickert]

Description

The docstring for v0.9.0 of tfp.stats.percentile states that

https://github.com/tensorflow/probability/blob/356cfddef026b3339b8f2a81e600acd2ff8e22b4/tensorflow_probability/python/stats/quantiles.py#L415-L416

From this description (specifically the "from the minimum to the maximum" part) I would assume that this means that for

>>> import tensorflow as tf
>>> a = tf.constant([[10.0, 7.0, 4.0], [3.0, 2.0, 1.0]])
>>> a
<tf.Tensor: shape=(2, 3), dtype=float32, numpy=
array([[10.,  7.,  4.],
       [ 3.,  2.,  1.]], dtype=float32)>

the sorted version of it that percentile evaluates on would be in _ascending_ order

>>> y = tf.reshape(a, [-1])
>>> tf.sort(y)
<tf.Tensor: shape=(6,), dtype=float32, numpy=array([ 1.,  2.,  3.,  4.,  7., 10.], dtype=float32)>

However, as I learned in this Stack Overflow question (Different results with NumPy percentile and TensorFlow percentile for “nearest” interpolation method)

the sort order used is actually _descending_

https://github.com/tensorflow/probability/blob/356cfddef026b3339b8f2a81e600acd2ff8e22b4/tensorflow_probability/python/stats/quantiles.py#L884

which gives

>>> sorted_y, _ = tf.math.top_k(y, k=tf.shape(y)[-1])
>>> sorted_y
<tf.Tensor: shape=(6,), dtype=float32, numpy=array([10.,  7.,  4.,  3.,  2.,  1.], dtype=float32)>

Problem

The above means that when the "nearest" interpolation method is used for tfp.stats.percentile this example gives

>>> import tensorflow as tf
>>> import tensorflow_probability as tfp
>>> a = tf.constant([[10.0, 7.0, 4.0], [3.0, 2.0, 1.0]])
>>> q = 50
>>> tfp.stats.percentile(a, q, interpolation="nearest")
<tf.Tensor: shape=(), dtype=float32, numpy=4.0>

where the docstring would make one think that the result would be

<tf.Tensor: shape=(), dtype=float32, numpy=3.0>

Additionally, as is outlined in my Stack Overflow question, this behavior gives different results than NumPy's numpy.percentile, which is confusing given that tfp.stats.percentile's docstring directly compares itself to numpy.percentile

https://github.com/tensorflow/probability/blob/356cfddef026b3339b8f2a81e600acd2ff8e22b4/tensorflow_probability/python/stats/quantiles.py#L428

Proposed Change

Given that the it seems that tfp.stats.percentile intentionally gives congruent behavior to numpy.percentile, it would be good to change the sort order used from _descending_ to _ascending_.

So maybe something like changing

https://github.com/tensorflow/probability/blob/356cfddef026b3339b8f2a81e600acd2ff8e22b4/tensorflow_probability/python/stats/quantiles.py#L884

to

sorted_, _ = tf.math.top_k(tensor, k=tf.shape(tensor)[-1])
sorted_ = tf.reverse(sorted_, [-1])

This would also make the code base more readable, as @langmore notes in the comments of _get_indices

https://github.com/tensorflow/probability/blob/356cfddef026b3339b8f2a81e600acd2ff8e22b4/tensorflow_probability/python/stats/quantiles.py#L550-L552

Additional Notes

This affects pyhf Issue 815 (I'm on the dev team).

cc @lukasheinrich @kratsg

[matthewfeickert]

If any of the tfp dev team can give thoughts on if this change would be something that would be considered it would be great to know. I'm a bit hesitant to start working on a PR if there was strong opposition to a change.

[csuter]

Thanks @matthewfeickert. This sounds like a bug/oversight on our part, but this behavior has been around for quite a while which means any existing uses of it might be depending on it and would be broken by a change. To change it safely, we'd probably need to introduce a flag that allows new users to opt into the "right" behavior, and then deprecate the old default over a period of 1-2 release cycles. @jvdillon @langmore any thoughts?

[matthewfeickert]

but this behavior has been around for quite a while which means any existing uses of it might be depending on it and would be broken by a change. To change it safely, we'd probably need to introduce a flag that allows new users to opt into the "right" behavior, and then deprecate the old default over a period of 1-2 release cycles.

Thanks very much for your thoughts and feedback, @csuter. Yup. I wasn't under any illusion that this would be a "quick" fix. It is good to hear though that you think this might be a possibility for change in the future.

I'll be watching this Issue if there is anything that the pyhf dev team can do to be helpful.

[csuter]

Oh cool, I didn't realize pyhf was using TFP. I see you're using it in the TF backend, mostly for [log_]probs. FYI, we are in the process of supporting using TFP with JAX or numpy as our backend. Most distributions are working now, certainly Normal and Poisson so you might be able to unify your sources of those bits, if that's of any interest :) Simple demo notebook here.

[matthewfeickert]

FYI, we are in the process of supporting using TFP with JAX or numpy as our backend.

That sounds pretty interesting. :) Thanks very much for this heads up!

[jvdillon]

I apologize for the bug. Thanks for the catch and nice explanation.

Im all for making the breaking change which does the "right" thing, ie np
behavior. If that's intolerable then I vote for Chris's proposal. (Eg,
"nearest_numpy".)

On Wed, Apr 15, 2020 at 3:48 PM Matthew Feickert [e-mail]
wrote:

FYI, we are in the process of supporting using TFP with JAX or numpy as
our backend.

That sounds pretty interesting. :) Thanks very much for this heads up!

—
You are receiving this because you were mentioned.
Reply to this email directly, view it on GitHub
https://github.com/tensorflow/probability/issues/864#issuecomment-614317721,
or unsubscribe
https://github.com/notifications/unsubscribe-auth/AAIVTNRKNZHCEOCN2NACQJTRMY2TBANCNFSM4L7PC74A
.

[matthewfeickert]

I'm all for making the breaking change which does the "right" thing, ie np behavior. If that's intolerable then I vote for Chris's proposal.

@csuter @jvdillon Do you have any further thoughts on what is the preferred path forward? If so, regardless of what they are, would it be helpful if I opened a PR? Or given other work that the TFP team is doing for the next release it is better to hold off and have this get picked up by a team member later on?

[csuter]

Sorry -- I misunderstood (overestimated) the extent of the change of behavior. Sounds like this is just fixing a broken thing, which we should just do. @matthewfeickert we'd welcome your generous offer of a PR! Thanks for the careful writeup and for diving in and diagnosing this subtle issue!

[matthewfeickert]

Sounds like this is just fixing a broken thing, which we should just do.

Fantastic! Thanks very much for this information.

we'd welcome your generous offer of a PR!

Great! :+1: I'll try to get to it before the end of the week if things don't get too hectic.