betaincinv skips the root finding when y == 0.5, so Beta medians are wrong

Abierto
  • ¿Algo no está bien?
#2,034 0 comentarios 0 reacciones 0 asignados Ver en GitHub

Nadie ha tomado este issue todavía.

Evaluación

Dificultad
3/5
Tiempo estimado
1-2 días
Aptitud para principiantes
68/100
Tipo de issue
Error
Claridad
Bien especificado
Estado de actividad
Activo
Stack tecnológico
python
Área
backend

Línea de trabajo

Start in tensorflow_probability/python/math/special.py at _betaincinv_computation and run the issue's reproduction to observe the incorrect median results. Add a regression test for y == 0.5, including the symmetric Beta case, and verify that the quantile results are correct without regressing other inputs.

Escrito por el modelo de indexación a partir del texto del issue.

Descripción

tfp.math.betaincinv(a, b, y) returns its initial approximation, without any root-finding iterations, when y == 0.5. So every median computed through it is off, including tfd.Beta(a, b).quantile(0.5). In float64 this happens for all a, b. In float32 it happens for part of the a < 1, b < 1 region, where the bracket is halved: (0.5, 0.5) gives 0.6185 and (0.9, 0.5) gives 0.7783, against 0.5 and 0.7210.

System info
  • tfp-nightly 0.26.0.dev20261008 with jax 0.10.2, and tensorflow-probability 0.25.0 with jax 0.4.30. Same output on both.
  • Linux (WSL2), CPU. The script downloads nothing.
Reproduction
import jax

jax.config.update("jax_enable_x64", True)
import jax.numpy as jnp
import scipy.stats as st
from tensorflow_probability.substrates import jax as tfp

for a, b in [(2.0, 3.0), (10.0, 1.0), (0.5, 0.5)]:
    got = tfp.math.betaincinv(jnp.float64(a), jnp.float64(b), jnp.float64(0.5))
    quantile = tfp.distributions.Beta(jnp.float64(a), jnp.float64(b)).quantile(jnp.float64(0.5))
    near = tfp.math.betaincinv(jnp.float64(a), jnp.float64(b), jnp.float64(0.5 - 1e-12))
    print(f"a={a} b={b}  betaincinv(y=0.5)={float(got):.6f}  Beta.quantile(0.5)={float(quantile):.6f}  "
          f"betaincinv(y=0.5-1e-12)={float(near):.6f}  scipy={st.beta(a, b).ppf(0.5):.6f}")

Current result:

a=2.0 b=3.0  betaincinv(y=0.5)=0.390087  Beta.quantile(0.5)=0.390087  betaincinv(y=0.5-1e-12)=0.385728  scipy=0.385728
a=10.0 b=1.0  betaincinv(y=0.5)=0.906306  Beta.quantile(0.5)=0.906306  betaincinv(y=0.5-1e-12)=0.933033  scipy=0.933033
a=0.5 b=0.5  betaincinv(y=0.5)=0.381514  Beta.quantile(0.5)=0.381514  betaincinv(y=0.5-1e-12)=0.499537  scipy=0.500000

(Beta(10, 1) has median 0.5 ** (1 / 10) = 0.933033.)

Root cause

In _betaincinv_computation the while loop starts with this stop mask (special.py#L753):

      loop_vars=(
          tf.equal(y, initial_low) | tf.equal(y, initial_high),
          initial_low,
          initial_high,
          initial_candidate),

initial_low and initial_high bracket the solution x, not the probability y. In float64 initial_high is always 0.5 (in float32 it is 0.5 in the small a, b branch), so y == 0.5 marks the element as finished before the first iteration and initial_candidate comes back as the result. The y == 0 and y == 1 cases are already handled after the loop, so I think the intended mask is:

          tf.equal(y, zero) | tf.equal(y, one),

With that one-line change the first two rows give 0.385728 and 0.933033, matching scipy. The mask only differs when y is 0.5, so other outputs are unchanged: 4000 random (a, b, y) triples came out bitwise identical before and after.

A second, smaller problem

For a == b the root at y = 0.5 sits right at the 0.5 bracket edge. Halley steps land just past high, the bisection fallback only halves the remaining gap, and after the 8 float64 iterations (0.5, 0.5) ends at 0.499537. This is the y=0.5-1e-12 column above, which goes through the full loop today. The one-line fix doesn't change this. Returning 0.5 directly when betainc(a, b, 0.5) is within tolerance of y would cover it, but I haven't tested that.

Downstream

NumPyro's Beta.icdf and StudentT.icdf call this function. StudentT.icdf passes 1 - |1 - 2q|, which is 0.5 at both quartiles, so in x64 dist.StudentT(3.0).icdf(0.75) gives 0.664628 instead of 0.764892. TFP's own StudentT.quantile isn't affected, since it goes through stdtrit.

I searched the issues for betaincinv, Beta quantile, StudentT quantile and quantile 0.5. The closest are #1581 and #1625 (2022, which added the current implementation), and neither mentions this. Happy to open a PR with the change and a regression test if that would help.

I used Claude Code to help trace this. I ran the reproduction and checked the numbers myself.

Lenguaje dominante
Jupyter Notebook
Estrellas
4.4k
Forks
1.1k
Sin PR fusionados en 30 d

Preparar el entorno

Primeros pasos

  1. Lee el issue completo y luego la guía de contribución del proyecto.
  2. Comenta en el issue que vas a ocuparte — evita que dos personas hagan lo mismo.
  3. Haz un fork del repositorio y trabaja en una rama.
  4. Abre un pull request que haga referencia al número del issue.

Más de tensorflow/probability

Todos los issues de tensorflow/probability

Issues similares

Más issues de Backend & API Design

Recibe los nuevos issues en tu correo

Un resumen breve de issues de GitHub para principiantes.