#34975 BUG confusion_matrix_at_thresholds / roc_auc_score lose precision on plain numpy when y_score is float32 (regression from #34817)

Describe the bug

confusion_matrix_at_thresholds (and everything built on it: roc_curve, precision_recall_curve, det_curve, roc_auc_score) silently loses integer precision on plain numpy — not just the restricted array-API devices that #34817 was fixing — whenever the caller's y_score array is float32. Model scores/predictions are very commonly stored as float32 (e.g. to save memory), so this is not an exotic input.

This looks like a regression introduced this week by #34817 (merged 2026-08-31), which switched the unweighted cumulative-sum path to compute exact int64 counts and then cast them to output_dtype:

https://github.com/scikit-learn/scikit-learn/blob/main/sklearn/metrics/_ranking.py#L1028-L1039

y_true_int = xp.astype(y_true, xp.int64)
tps_int = xp.cumulative_sum(y_true_int, dtype=xp.int64)[threshold_idxs]
fps_int = (xp.astype(threshold_idxs, xp.int64) + 1) - tps_int

output_dtype = (
    y_score.dtype
    if hasattr(y_score, "dtype") and xp.isdtype(y_score.dtype, "real floating")
    else _max_precision_float_dtype(xp, device)
)
tps = xp.astype(tps_int, output_dtype)
fps = xp.astype(fps_int, output_dtype)

output_dtype is set to y_score.dtype whenever y_score is real-floating, instead of the device's max-precision float dtype (_max_precision_float_dtype(xp, device), which is always float64 on plain numpy). So on ordinary numpy — where float64 is always available — if y_score happens to be float32, the exact int64 counts get needlessly downcast to float32, which cannot represent integers exactly above 2**24 (~16.7M). The weighted branch right below doesn't have this problem: it always uses max_float_dtype regardless of y_score's dtype.

This is exactly the precision-loss bug class #34817/#34813 set out to fix, just reintroduced on the ordinary numpy path whenever the input y_score happens to be float32 — a case the fix's own test suite doesn't cover (the added test test_confusion_matrix_at_thresholds_float32_only_unweighted_large_n only exercises the array_api_strict "no float64" restricted-device scenario, not plain numpy with a float32 y_score).

Steps/Code to Reproduce

import numpy as np
from sklearn.metrics import confusion_matrix_at_thresholds, roc_auc_score

rng = np.random.RandomState(0)
n = 20_000_000
y_true = (rng.rand(n) < 0.9).astype(np.int32)
y_score64 = rng.rand(n).astype(np.float64)
y_score32 = y_score64.astype(np.float32)  # e.g. model scores stored as float32

_, _, _, tps64, _ = confusion_matrix_at_thresholds(y_true, y_score64)
_, _, _, tps32, _ = confusion_matrix_at_thresholds(y_true, y_score32)

print("tps64 dtype:", tps64.dtype, "last value:", tps64[-1])
print("tps32 dtype:", tps32.dtype, "last value:", tps32[-1])
print("n_pos (exact):", y_true.sum())

print("roc_auc_score (float64 y_score):", roc_auc_score(y_true, y_score64))
print("roc_auc_score (float32 y_score):", roc_auc_score(y_true, y_score32))

Expected Results

tps should be exact regardless of y_score's dtype, on a device where float64 is available (matching the weighted branch's behavior, and matching pre-#34817 behavior):

tps64 dtype: float64 last value: 18000923.0
tps32 dtype: float64 last value: 18000923.0
n_pos (exact): 18000923

Actual Results

tps64 dtype: float64 last value: 18000923.0
tps32 dtype: float32 last value: 18000924.0
n_pos (exact): 18000923
roc_auc_score (float64 y_score): 0.500012352595143
roc_auc_score (float32 y_score): 0.5000123248180354

The float32-y_score run silently gives a wrong (off-by-one, and for larger n increasingly wrong) true-positive count and a slightly different roc_auc_score, purely because of y_score's dtype — with no warning and no restricted-device context involved.

Versions

Reproduced on current main (post-#34817, i.e. including the fix for #34813/#34827) with plain numpy — no array_api_dispatch context, no array-API device involved.

I dug into this while looking for an unrelated fix and haven't opened a PR — happy to put one together (the obvious fix is to use _max_precision_float_dtype(xp, device) unconditionally for output_dtype, mirroring the weighted branch, since it's already computed as the fallback) if that direction looks right to a maintainer.