chisquared entropy/logpdf lose float64 precision — pt.log(2) folds at float32

Summary

chisquared.entropy and chisquared.logpdf use pt.log(2) with a Python int literal. PyTensor stores 2 as int8, and log upgrades int8/int16 inputs to float32, so the constant folds to 0.6931471824645996 instead of the float64 log(2) = 0.6931471805599453. That ~1.9e-9 error leaks into both functions.

Affected code

• chisquared.py:54 (entropy): ... + pt.log(2) + ...
• chisquared.py:90 (logpdf): ... - (nu * pt.log(2)) / 2

Reproduction

import numpy as np, pytensor.tensor as pt
from scipy import stats
from pytensor_distributions import chisquared as C

print(float(pt.log(2).eval()))          # 0.6931471824645996  (float64 log(2) = 0.6931471805599453)

for nu in [2.0, 10.0, 100.0]:
    e  = float(C.entropy(pt.constant(nu, "float64")).eval())
    lp = float(C.logpdf(pt.constant(nu, "float64"), pt.constant(nu, "float64")).eval())
    print(nu, abs(e - stats.chi2(nu).entropy()), abs(lp - stats.chi2(nu).logpdf(nu)))

Measured error vs scipy:

entropy is a flat ~1.9e-9; logpdf's scales with nu (the nu·log(2) term), so it grows with the degrees of freedom and accumulates across a summed log-likelihood.

Root cause (upstream, not really a chisquared bug)

This is the "aggressive downcast of literals" root cause (pymc-devs/pytensor#1073) hitting direct user code: pt.log(<int>) → int8 → float32. A maintainer has confirmed #1073 is the underlying issue (comment):

the underlying issue is the same old aggressive downcast of literals: #1073

The obvious pt.log(2.0) does not work either: the Python float 2.0 autocasts to float32, so it's still 0.6931471824645996. Only a NumPy float64 constant is correct.

The real fix belongs upstream (#1073) — once PyTensor stops downcasting Python literals, pt.log(2) will be exact and nothing here needs changing. That root fix (pymc-devs/pytensor#1917) is still an open draft, so a stopgap may be warranted in the meantime.

Workaround (stopgap only — not the real fix)

Until #1073 lands, sidestep the downcast by feeding log a NumPy float64 constant, matching the pattern already used in this repo (genpareto.median, and the merged asymmetriclaplace.median fix in #60):

import numpy as np
# entropy:  ... + np.log(2.0) + ...
# logpdf:   ... - (nu * np.log(2.0)) / 2

This is a local workaround, not a fix for the underlying downcast; it should be revisited (and ideally reverted) once #1073 is resolved.

Scope

The same pt.log(<int>) pattern appears in ~9 modules (inversegamma, moyal, wishart, halfcauchy, weibull, halfstudentt, skew_studentt, polyagamma, plus chisquared) — all the same upstream symptom, all sidesteppable the same way pending #1073. Worth a repo-wide sweep if a stopgap is desired.

Start in chisquared.py at the entropy line 54 and logpdf line 90, then inspect the existing genpareto.median pattern and the asymmetriclaplace.median fix referenced in the issue. Run the provided reproduction against scipy; done means both functions avoid the measured float32 precision error while the upstream PyTensor issue remains unresolved.