From dc60819275f6b8d2bd75becfc6df5accb2cede88 Mon Sep 17 00:00:00 2001 From: Seongho Bae Date: Mon, 27 Jul 2026 08:04:16 +0900 Subject: [PATCH 1/2] Add predict_rating and predict_rating_multi (PlayerRatings predict.rating) Rust core predict_rating_two/predict_rating_multi in mlsirm-core scaling (Elo logistic, deviation-shrunk Glicko-family with qip3 = 3(ln10/400/pi)^2, EloM rowmean branch with optional min-tie placing), PyO3 bindings, Python wrappers, exact-oracle anchor tests P1-P9, 7-mutant EXECUTED kill map, MC-500 invariants, TestPredict pytest coverage, CHANGELOG entry. Normative source: CRAN PlayerRatings 1.1-0 R/ratings.R lines 1056-1133 (READ). REDUCED-SCOPE: index-based players, per-game/scalar gamma. Co-authored-by: Copilot App <223556219+Copilot@users.noreply.github.com> --- CHANGELOG.md | 1 + crates/fast-mlsirm-py/src/lib.rs | 73 ++++ crates/mlsirm-core/src/scaling.rs | 291 +++++++++++++ python/fast_mlsirm/__init__.py | 4 + python/fast_mlsirm/scaling.py | 232 +++++++++++ tests/test_paper_features.py | 134 ++++++ tests/unit/scaling_tests.rs | 663 ++++++++++++++++++++++++++++++ 7 files changed, 1398 insertions(+) diff --git a/CHANGELOG.md b/CHANGELOG.md index bc9855d5a..9c6c685e0 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -119,6 +119,7 @@ ### Added +- `predict_rating` / `predict_rating_multi` game-outcome prediction from fitted ratings (CRAN PlayerRatings 1.1-0 `predict.rating` `R/ratings.R` lines 1056–1133 — READ and normative; no journal paper exists for this dispatch, CRAN package provenance only): Elo logistic branch, deviation-shrunk Glicko/Glicko-2/Stephenson branch (`qip3 = 3(ln10/400π)²`, joint shrink over BOTH players' squared deviations), and multi-player EloM branch (`(rating − rowmean)/40`, na.rm row means, optional min-tie placing ranks with NaN kept). R semantics preserved: strict `games < tng` unrated cutoff, `trat` replacement of ALL missing extracted values (unmatched, low-games, stored-NA), `pred >= thresh` binarization with NaN propagation. REDUCED-SCOPE vs R: index-based (−1 = unmatched; caller does name matching), per-game/scalar gamma only. Exact-oracle fixtures P1–P9 and a 7-mutant EXECUTED kill map (incl. both branches' tng comparisons). - `fide_rating` FIDE-style Elo ratings (CRAN PlayerRatings 1.1-0 `fide()` `R/ratings.R` lines 125–272 + `kfide()` lines 959–972 — READ and normative; no journal paper exists for this variant, CRAN package provenance only): per-period batch Elo with the kfide K-factor schedule (K = kv[0] elite / kv[1] ≥30 games / kv[2] novice, evaluated from PERIOD-START state), sticky elite flag set from POST-update ratings ≥ 2400, and per-player running mean of POST-update opponent ratings. REDUCED-SCOPE vs R: no status/history frames, kfide-only K schedule, self-play rejected, thresholds 30/2400 hard-coded. kv=(k,k,k) reduces bitwise to `elo_rating(kfac=k)` (MC-500 anchor); exact-oracle fixtures F1–F5 and a 5-mutant EXECUTED kill map. - `metrics_rating` prediction-quality metrics for binary-outcome forecasts (binomial deviance on capped predictions, RMSE/MAE on raw predictions, optional 0.5-baseline scaling), a Rust reimplementation of CRAN PlayerRatings 1.1-0 `metrics()` with its cap quirk and elementwise NaN semantics preserved; exact-oracle anchor tests, mutation-kill map, and MC-500 invariants. - **Multiplayer Elo rating (CRAN PlayerRatings 1.1-0 `elom()` diff --git a/crates/fast-mlsirm-py/src/lib.rs b/crates/fast-mlsirm-py/src/lib.rs index 7d2ddd264..b866b2d50 100644 --- a/crates/fast-mlsirm-py/src/lib.rs +++ b/crates/fast-mlsirm-py/src/lib.rs @@ -5741,6 +5741,77 @@ fn fide_rating( Ok(d.into()) } +/// Predicted game outcomes from fitted ratings, PlayerRatings +/// `predict.rating` two-player branches (see +/// `mlsirm_core::scaling::predict_rating_two`). `white`/`black` are player +/// indices with -1 = unmatched. Pass `deviations` for the +/// Glicko/Glicko-2/Stephenson deviation-shrunk branch. +#[pyfunction] +#[pyo3(signature = (ratings, deviations, games, white, black, gamma, tng, trat_rating, trat_deviation, thresh))] +#[allow(clippy::too_many_arguments)] +fn predict_rating_two( + py: Python<'_>, + ratings: PyReadonlyArray1<'_, f64>, + deviations: Option>, + games: PyReadonlyArray1<'_, u64>, + white: PyReadonlyArray1<'_, i64>, + black: PyReadonlyArray1<'_, i64>, + gamma: PyReadonlyArray1<'_, f64>, + tng: u64, + trat_rating: Option, + trat_deviation: Option, + thresh: Option, +) -> PyResult>> { + let dev_slice = deviations.as_ref().map(|d| d.as_slice()).transpose()?; + let trat = trat_rating.map(|t1| (t1, trat_deviation.unwrap_or(f64::NAN))); + let out = mlsirm_core::scaling::predict_rating_two( + ratings.as_slice()?, + dev_slice, + games.as_slice()?, + white.as_slice()?, + black.as_slice()?, + gamma.as_slice()?, + tng, + trat, + thresh, + ) + .map_err(PyValueError::new_err)?; + Ok(PyArray1::from_slice(py, &out).into()) +} + +/// Predicted expected scores for multi-player (EloM) events, PlayerRatings +/// `predict.rating` EloM branch (see +/// `mlsirm_core::scaling::predict_rating_multi`). `players` is flattened +/// row-major nr x np with -1 = empty seat; `placing` returns min-tie ranks. +#[pyfunction] +#[pyo3(signature = (ratings, games, players, nr, np, tng, trat, placing))] +#[allow(clippy::too_many_arguments)] +fn predict_rating_multi( + py: Python<'_>, + ratings: PyReadonlyArray1<'_, f64>, + games: PyReadonlyArray1<'_, u64>, + players: PyReadonlyArray1<'_, i64>, + nr: usize, + np: usize, + tng: u64, + trat: Option, + placing: bool, +) -> PyResult>> { + let out = mlsirm_core::scaling::predict_rating_multi( + ratings.as_slice()?, + games.as_slice()?, + players.as_slice()?, + nr, + np, + tng, + trat, + placing, + ) + .map_err(PyValueError::new_err)?; + Ok(PyArray1::from_slice(py, &out).into()) +} + + /// GPCM/nominal softmax cell log-probabilities at one node (parity surface for /// the NumPy `category_logprobs` reference). #[pyfunction] @@ -8095,6 +8166,8 @@ fn fast_mlsirm_core(m: &Bound<'_, PyModule>) -> PyResult<()> { m.add_function(wrap_pyfunction!(elom_rating, m)?)?; m.add_function(wrap_pyfunction!(metrics_rating, m)?)?; m.add_function(wrap_pyfunction!(fide_rating, m)?)?; + m.add_function(wrap_pyfunction!(predict_rating_two, m)?)?; + m.add_function(wrap_pyfunction!(predict_rating_multi, m)?)?; m.add_function(wrap_pyfunction!(circle_arc_middle_anchor, m)?)?; m.add_function(wrap_pyfunction!(loglinear_smooth, m)?)?; m.add_function(wrap_pyfunction!(person_fit_stat, m)?)?; diff --git a/crates/mlsirm-core/src/scaling.rs b/crates/mlsirm-core/src/scaling.rs index 8a88c80e4..1d918dc65 100644 --- a/crates/mlsirm-core/src/scaling.rs +++ b/crates/mlsirm-core/src/scaling.rs @@ -3701,6 +3701,297 @@ pub fn fide_rating( }) } +/// Game-outcome prediction from fitted ratings — two-player branches of +/// `predict.rating` from CRAN PlayerRatings 1.1-0 (`R/ratings.R` lines +/// 1056-1133; source READ). No journal paper exists for this dispatch +/// function; the CRAN R source is the normative reference (provenance as +/// for `stephenson_rating`). Two branches are implemented here: +/// +/// * Elo branch (fitted `type == "Elo"`, produced by `elo()`/`fide()`): +/// `pred = 1 / (1 + 10^((brat - wrat - gamma)/400))` (R line 1118). +/// * Deviation branch (`Glicko`/`Glicko-2`/`Stephenson`), engaged when +/// `deviations` is supplied: with `qv = ln(10)/400` and +/// `qip3 = 3 (qv/pi)^2`, +/// `vec = 1/sqrt(1 + qip3 (wdev^2 + bdev^2))` and +/// `pred = 1 / (1 + 10^(vec (brat - wrat - gamma)/400))` +/// (R lines 1071, 1119-1121; note the joint shrink sums BOTH players' +/// squared deviations). +/// +/// Pre-processing (verified against R 1069-1114 and the executed oracle): +/// players with `games < tng` (strict `<`) have their stored rating (and +/// deviation) set to NA before extraction; when `trat` is supplied it +/// then replaces ALL extracted-NA values — unmatched players (index +/// sentinel `-1`), low-games players, and matched players whose stored +/// rating/deviation is already NaN. Without `trat`, NA propagates to a +/// NaN prediction. `thresh` (R 1127-1128) maps `pred >= thresh` to 1 else +/// 0, leaving NaN predictions NaN. +/// +/// REDUCED SCOPE relative to R: index-based (caller does name matching; +/// `-1` = unmatched), per-game `gamma` (scalar recycling is a wrapper +/// concern), no data-frame interface or `object$type` dispatch. +/// +/// All formulas were verified against an executed exact/float oracle +/// transcribed from the R source; claims are limited to that source. +pub fn predict_rating_two( + ratings: &[f64], + deviations: Option<&[f64]>, + games: &[u64], + white: &[i64], + black: &[i64], + gamma: &[f64], + tng: u64, + trat: Option<(f64, f64)>, + thresh: Option, +) -> Result, String> { + let n = ratings.len(); + if n < 2 || n > 10_000 { + return Err(format!( + "predict_rating_two: number of players must be in 2..=10000, got {n}" + )); + } + if games.len() != n { + return Err(format!( + "predict_rating_two: games length {} != number of players {n}", + games.len() + )); + } + if let Some(dev) = deviations { + if dev.len() != n { + return Err(format!( + "predict_rating_two: deviations length {} != number of players {n}", + dev.len() + )); + } + if dev.iter().any(|d| d.is_infinite()) { + return Err( + "predict_rating_two: deviations must not be infinite (NaN = missing)".to_string(), + ); + } + } + if ratings.iter().any(|r| r.is_infinite()) { + return Err("predict_rating_two: ratings must not be infinite (NaN = missing)".to_string()); + } + let g = white.len(); + if g == 0 { + return Err("predict_rating_two: at least one game row is required".to_string()); + } + if black.len() != g || gamma.len() != g { + return Err(format!( + "predict_rating_two: white/black/gamma lengths must match, got {}/{}/{}", + g, + black.len(), + gamma.len() + )); + } + if gamma.iter().any(|x| !x.is_finite()) { + return Err("predict_rating_two: gamma must be finite".to_string()); + } + if let Some((t1, t2)) = trat { + if !t1.is_finite() { + return Err("predict_rating_two: trat rating must be finite".to_string()); + } + if deviations.is_some() && !t2.is_finite() { + return Err("predict_rating_two: trat deviation must be finite".to_string()); + } + } + if let Some(t) = thresh { + if !t.is_finite() { + return Err("predict_rating_two: thresh must be finite".to_string()); + } + } + for (w, b) in white.iter().zip(black.iter()) { + for idx in [*w, *b] { + if idx < -1 || idx >= n as i64 { + return Err(format!( + "predict_rating_two: player index {idx} out of range (-1 = unmatched, else 0..{n})" + )); + } + } + if *w >= 0 && w == b { + return Err("predict_rating_two: self-play rows are not allowed".to_string()); + } + } + // R 1088-1093: low-games rows lose their stored values before extraction. + let effective = |p: i64, table: &[f64]| -> f64 { + if p < 0 { + f64::NAN + } else if games[p as usize] < tng { + f64::NAN + } else { + table[p as usize] + } + }; + let qv = std::f64::consts::LN_10 / 400.0; + let qip3 = 3.0 * (qv / std::f64::consts::PI) * (qv / std::f64::consts::PI); + let mut preds = Vec::with_capacity(g); + for k in 0..g { + let mut wrat = effective(white[k], ratings); + let mut brat = effective(black[k], ratings); + if let Some((t1, _)) = trat { + // R 1097-1098: trat replaces ALL extracted NA ratings. + if wrat.is_nan() { + wrat = t1; + } + if brat.is_nan() { + brat = t1; + } + } + let pred = if let Some(dev) = deviations { + let mut wdev = effective(white[k], dev); + let mut bdev = effective(black[k], dev); + if let Some((_, t2)) = trat { + if wdev.is_nan() { + wdev = t2; + } + if bdev.is_nan() { + bdev = t2; + } + } + if wrat.is_nan() || brat.is_nan() || wdev.is_nan() || bdev.is_nan() { + f64::NAN + } else { + let vec = 1.0 / (1.0 + qip3 * (wdev * wdev + bdev * bdev)).sqrt(); + 1.0 / (1.0 + 10f64.powf(vec * (brat - wrat - gamma[k]) / 400.0)) + } + } else if wrat.is_nan() || brat.is_nan() { + f64::NAN + } else { + 1.0 / (1.0 + 10f64.powf((brat - wrat - gamma[k]) / 400.0)) + }; + preds.push(pred); + } + if let Some(t) = thresh { + // R 1127-1128: as.numeric(preds >= thresh); NA stays NA. + for p in preds.iter_mut() { + if !p.is_nan() { + *p = if *p >= t { 1.0 } else { 0.0 }; + } + } + } + Ok(preds) +} + +/// Multi-player (EloM) branch of `predict.rating` from CRAN PlayerRatings +/// 1.1-0 (`R/ratings.R` lines 1103-1105, 1123-1125, 1129-1130; source +/// READ; provenance as for `predict_rating_two`). +/// +/// Per event row: `pred = (rat - rowmean)/40` where `rowmean` is the mean +/// over the row's non-NaN seat ratings (`rowMeans(rats, na.rm=TRUE)`, +/// all-NaN row -> NaN row). Pre-processing (tng/trat) as in +/// `predict_rating_two`. With `placing = true` the per-row predictions +/// are replaced by `rank(-preds, na.last="keep", ties.method="min")`: +/// rank 1 = highest prediction, ties share the minimum rank, NaN seats +/// keep NaN. +/// +/// `players` is row-major `nr x np` with `-1` = empty/unmatched seat. +/// All formulas were verified against an executed oracle transcribed +/// from the R source. +pub fn predict_rating_multi( + ratings: &[f64], + games: &[u64], + players: &[i64], + nr: usize, + np: usize, + tng: u64, + trat: Option, + placing: bool, +) -> Result, String> { + let n = ratings.len(); + if n < 2 || n > 10_000 { + return Err(format!( + "predict_rating_multi: number of players must be in 2..=10000, got {n}" + )); + } + if games.len() != n { + return Err(format!( + "predict_rating_multi: games length {} != number of players {n}", + games.len() + )); + } + if ratings.iter().any(|r| r.is_infinite()) { + return Err( + "predict_rating_multi: ratings must not be infinite (NaN = missing)".to_string(), + ); + } + if nr == 0 { + return Err("predict_rating_multi: at least one event row is required".to_string()); + } + if np < 2 || np > 1000 { + return Err(format!( + "predict_rating_multi: seats per event must be in 2..=1000, got {np}" + )); + } + if players.len() != nr * np { + return Err(format!( + "predict_rating_multi: players length {} != nr*np = {}", + players.len(), + nr * np + )); + } + if let Some(t) = trat { + if !t.is_finite() { + return Err("predict_rating_multi: trat must be finite".to_string()); + } + } + for idx in players { + if *idx < -1 || *idx >= n as i64 { + return Err(format!( + "predict_rating_multi: player index {idx} out of range (-1 = empty, else 0..{n})" + )); + } + } + let mut out = vec![f64::NAN; nr * np]; + let mut row_rat = vec![f64::NAN; np]; + for i in 0..nr { + for s in 0..np { + let p = players[i * np + s]; + let mut v = if p < 0 || games[p as usize] < tng { + f64::NAN + } else { + ratings[p as usize] + }; + if let Some(t) = trat { + if v.is_nan() { + v = t; + } + } + row_rat[s] = v; + } + // rowMeans(rats, na.rm = TRUE) — R 1123. + let mut sum = 0.0; + let mut cnt = 0usize; + for v in row_rat.iter() { + if !v.is_nan() { + sum += v; + cnt += 1; + } + } + let mean = if cnt > 0 { sum / cnt as f64 } else { f64::NAN }; + for s in 0..np { + if !row_rat[s].is_nan() { + out[i * np + s] = (row_rat[s] - mean) / 40.0; + } + } + if placing { + // rank(-preds, na.last="keep", ties.method="min") — R 1129-1130: + // rank of v = 1 + count of non-NaN row values strictly greater. + for s in 0..np { + row_rat[s] = out[i * np + s]; + } + for s in 0..np { + if !row_rat[s].is_nan() { + let greater = row_rat + .iter() + .filter(|u| !u.is_nan() && **u > row_rat[s]) + .count(); + out[i * np + s] = 1.0 + greater as f64; + } + } + } + } + Ok(out) +} + #[cfg(test)] #[path = "../../../tests/unit/scaling_tests.rs"] mod tests; diff --git a/python/fast_mlsirm/__init__.py b/python/fast_mlsirm/__init__.py index 7bc8786df..67ded7a0c 100644 --- a/python/fast_mlsirm/__init__.py +++ b/python/fast_mlsirm/__init__.py @@ -41,6 +41,8 @@ metrics_rating as metrics_rating, fide_rating as fide_rating, FideResult as FideResult, + predict_rating as predict_rating, + predict_rating_multi as predict_rating_multi, ) from .cdm import fit_cdm as fit_cdm, CdmFit as CdmFit, fit_gdina as fit_gdina, GdinaFit as GdinaFit, validate_q_matrix as validate_q_matrix, QMatrixValidation as QMatrixValidation, gdina_wald_selection as gdina_wald_selection, WaldModelSelection as WaldModelSelection, fit_ho_cdm as fit_ho_cdm, HoCdmFit as HoCdmFit, fit_ho_gdina as fit_ho_gdina, HoGdinaFit as HoGdinaFit, fit_seq_gdina as fit_seq_gdina, SeqGdinaFit as SeqGdinaFit, fit_seq_gdina_qr as fit_seq_gdina_qr, SeqGdinaQrFit as SeqGdinaQrFit from .mixture import fit_mixture as fit_mixture, MixtureFit as MixtureFit @@ -274,6 +276,8 @@ "metrics_rating", "fide_rating", "FideResult", + "predict_rating", + "predict_rating_multi", "circle_arc_middle_anchor", "CircleArcResult", "fit_response_times", diff --git a/python/fast_mlsirm/scaling.py b/python/fast_mlsirm/scaling.py index efaded88a..3089faf93 100644 --- a/python/fast_mlsirm/scaling.py +++ b/python/fast_mlsirm/scaling.py @@ -1760,3 +1760,235 @@ def fide_rating(games, n_players, init=2200.0, kv=(10.0, 15.0, 30.0), gamma=None elite=np.asarray(res["elite"]), opponent=np.asarray(res["opponent"]), ) + + +def _predict_int_index_array(x, name, fname): + """Validate a player-index array: integers >= -1 (-1 = unmatched).""" + if isinstance(x, np.ma.MaskedArray): + raise ValueError(f"{fname}: masked arrays are not supported for {name}") + raw = x if isinstance(x, np.ndarray) else np.asarray(x, dtype=object) + if np.iscomplexobj(raw) or raw.dtype.kind == "b": + raise ValueError(f"{fname}: {name} must be integer indices, not complex/bool") + if raw.dtype == object: + for v in np.ravel(raw): + if v is None or isinstance( + v, (bool, np.bool_, str, bytes, np.datetime64, np.timedelta64) + ): + raise ValueError(f"{fname}: {name} contains a non-numeric value") + try: + arr = np.asarray(x, dtype=float) + except (TypeError, ValueError) as exc: + raise ValueError(f"{fname}: {name} is not numeric: {exc}") from None + if arr.ndim != 1: + raise ValueError(f"{fname}: {name} must be one-dimensional") + if not np.all(np.isfinite(arr)): + raise ValueError(f"{fname}: {name} must be finite") + if not np.all(arr == np.floor(arr)): + raise ValueError(f"{fname}: {name} must contain integers") + if np.any(arr < -1): + raise ValueError(f"{fname}: {name} indices must be >= -1 (-1 = unmatched)") + return arr.astype(np.int64) + + +def _predict_float_array(x, name, fname, allow_nan): + """Validate a float array; NaN optionally allowed (R NA), Inf rejected.""" + if isinstance(x, np.ma.MaskedArray): + raise ValueError(f"{fname}: masked arrays are not supported for {name}") + raw = x if isinstance(x, np.ndarray) else np.asarray(x, dtype=object) + if np.iscomplexobj(raw) or raw.dtype.kind == "b": + raise ValueError(f"{fname}: {name} must be real numeric, not complex/bool") + if raw.dtype == object: + for v in np.ravel(raw): + if isinstance( + v, (bool, np.bool_, str, bytes, np.datetime64, np.timedelta64) + ) or v is None: + raise ValueError(f"{fname}: {name} contains a non-numeric value") + try: + arr = np.asarray(x, dtype=float) + except (TypeError, ValueError) as exc: + raise ValueError(f"{fname}: {name} is not numeric: {exc}") from None + if arr.ndim != 1: + raise ValueError(f"{fname}: {name} must be one-dimensional") + if np.any(np.isinf(arr)): + raise ValueError(f"{fname}: {name} must not contain infinities") + if not allow_nan and np.any(np.isnan(arr)): + raise ValueError(f"{fname}: {name} must not contain NaN") + return arr + + +def _predict_scalar(x, name, fname): + """Validate a finite real scalar parameter.""" + import math + if isinstance(x, (bool, np.bool_)): + raise ValueError(f"{fname}: {name} must be real numeric, not bool") + raw = np.asarray(x) + if np.iscomplexobj(raw) or raw.dtype.kind == "b": + raise ValueError(f"{fname}: {name} must be real numeric, not complex/bool") + if raw.dtype == object and raw.ndim == 0 and isinstance( + raw.item(), (bool, np.bool_) + ): + raise ValueError(f"{fname}: {name} must be real numeric, not bool") + try: + v = float(x) + except (TypeError, ValueError) as exc: + raise ValueError(f"{fname}: {name} is not numeric: {exc}") from None + if not math.isfinite(v): + raise ValueError(f"{fname}: {name} must be finite") + return v + + +def predict_rating( + ratings, + games, + white, + black, + deviations=None, + gamma=30.0, + tng=15, + trat=None, + thresh=None, +): + """Predicted win probabilities for two-player games from fitted ratings. + + Implements the two-player branches of ``predict.rating`` from CRAN + PlayerRatings 1.1-0 (R/ratings.R lines 1056-1133; source read). Without + ``deviations`` this is the Elo branch (used for ``elo()``/``fide()`` + fits); with ``deviations`` it is the deviation-shrunk branch shared by + Glicko, Glicko-2, and Stephenson fits. Players are addressed by index; + ``-1`` marks an unmatched player (R's ``match()`` NA). Players with + ``games < tng`` (strict) are treated as unrated; ``trat`` (a scalar, or + a ``(rating, deviation)`` pair when ``deviations`` is supplied) + replaces ALL missing extracted values, matching R. ``thresh`` maps + predictions to 1 when ``pred >= thresh`` else 0 (NaN kept). ``gamma`` + is the per-game (or scalar, broadcast) first-player advantage. + """ + import math + + from .fitstats import _core_module + + fname = "predict_rating" + rat = _predict_float_array(ratings, "ratings", fname, allow_nan=True) + n = rat.shape[0] + if n < 2 or n > 10000: + raise ValueError(f"{fname}: number of players must be in 2..=10000, got {n}") + g_arr = _predict_float_array(games, "games", fname, allow_nan=False) + if np.any(g_arr < 0) or not np.all(g_arr == np.floor(g_arr)): + raise ValueError(f"{fname}: games must be nonnegative integers") + games_u64 = g_arr.astype(np.uint64) + w = _predict_int_index_array(white, "white", fname) + b = _predict_int_index_array(black, "black", fname) + dev = None + if deviations is not None: + dev = _predict_float_array(deviations, "deviations", fname, allow_nan=True) + ng = w.shape[0] + gam_raw = np.asarray(gamma) + if np.iscomplexobj(gam_raw) or gam_raw.dtype.kind == "b": + raise ValueError(f"{fname}: gamma must be real numeric, not complex/bool") + if gam_raw.ndim == 0: + gam = np.full(ng, _predict_scalar(gamma, "gamma", fname)) + else: + gam = _predict_float_array(gamma, "gamma", fname, allow_nan=False) + if gam.shape != (ng,): + raise ValueError( + f"{fname}: gamma must be a scalar or length-{ng} array, got {gam.shape}" + ) + tng_v = _predict_scalar(tng, "tng", fname) + if tng_v < 0 or tng_v != math.floor(tng_v): + raise ValueError(f"{fname}: tng must be a nonnegative integer") + trat_rating = None + trat_deviation = None + if trat is not None: + if dev is not None: + if not (isinstance(trat, (tuple, list)) and len(trat) == 2): + raise ValueError( + f"{fname}: trat must be a (rating, deviation) pair when " + "deviations are supplied" + ) + trat_rating = _predict_scalar(trat[0], "trat rating", fname) + trat_deviation = _predict_scalar(trat[1], "trat deviation", fname) + else: + if isinstance(trat, (tuple, list)): + if len(trat) != 1: + raise ValueError( + f"{fname}: trat must be a scalar (length 1) without deviations" + ) + trat = trat[0] + trat_rating = _predict_scalar(trat, "trat", fname) + thresh_v = None if thresh is None else _predict_scalar(thresh, "thresh", fname) + out = _core_module().predict_rating_two( + np.ascontiguousarray(rat), + None if dev is None else np.ascontiguousarray(dev), + np.ascontiguousarray(games_u64), + np.ascontiguousarray(w), + np.ascontiguousarray(b), + np.ascontiguousarray(gam), + int(tng_v), + trat_rating, + trat_deviation, + thresh_v, + ) + return np.asarray(out) + + +def predict_rating_multi( + ratings, + games, + players, + tng=15, + trat=None, + placing=False, +): + """Predicted expected scores (or placings) for multi-player EloM events. + + Implements the EloM branch of ``predict.rating`` from CRAN + PlayerRatings 1.1-0 (R/ratings.R lines 1103-1105, 1123-1125, + 1129-1130; source read): per event row, + ``pred = (rating - rowmean) / 40`` with the row mean over non-missing + seats (``rowMeans(rats, na.rm=TRUE)``). ``players`` is an ``(nr, np)`` + index matrix with ``-1`` = empty/unmatched seat. Players with + ``games < tng`` are treated as unrated; scalar ``trat`` replaces all + missing extracted ratings. With ``placing=True`` each row is replaced + by min-tie ranks of the predictions (rank 1 = highest; NaN kept), + matching R's ``rank(-preds, na.last="keep", ties.method="min")``. + """ + import math + + from .fitstats import _core_module + + fname = "predict_rating_multi" + rat = _predict_float_array(ratings, "ratings", fname, allow_nan=True) + n = rat.shape[0] + if n < 2 or n > 10000: + raise ValueError(f"{fname}: number of players must be in 2..=10000, got {n}") + g_arr = _predict_float_array(games, "games", fname, allow_nan=False) + if np.any(g_arr < 0) or not np.all(g_arr == np.floor(g_arr)): + raise ValueError(f"{fname}: games must be nonnegative integers") + games_u64 = g_arr.astype(np.uint64) + if isinstance(players, np.ma.MaskedArray): + raise ValueError(f"{fname}: masked arrays are not supported for players") + p_raw = players if isinstance(players, np.ndarray) else np.asarray(players, dtype=object) + if p_raw.ndim != 2: + raise ValueError(f"{fname}: players must be a 2-D (events, seats) matrix") + nr, np_seats = p_raw.shape + flat = _predict_int_index_array(np.ravel(p_raw), "players", fname) + if not (2 <= np_seats <= 1000): + raise ValueError( + f"{fname}: seats per event must be in 2..=1000, got {np_seats}" + ) + tng_v = _predict_scalar(tng, "tng", fname) + if tng_v < 0 or tng_v != math.floor(tng_v): + raise ValueError(f"{fname}: tng must be a nonnegative integer") + trat_v = None if trat is None else _predict_scalar(trat, "trat", fname) + if not isinstance(placing, (bool, np.bool_)): + raise ValueError(f"{fname}: placing must be a bool") + out = _core_module().predict_rating_multi( + np.ascontiguousarray(rat), + np.ascontiguousarray(games_u64), + np.ascontiguousarray(flat), + int(nr), + int(np_seats), + int(tng_v), + trat_v, + bool(placing), + ) + return np.asarray(out).reshape(nr, np_seats) diff --git a/tests/test_paper_features.py b/tests/test_paper_features.py index e45cdc783..8ea26dfae 100644 --- a/tests/test_paper_features.py +++ b/tests/test_paper_features.py @@ -10067,3 +10067,137 @@ def test_elo_reduction(self): assert f.ratings.tolist() == e.ratings.tolist() assert f.games.tolist() == e.games.tolist() assert f.lag.tolist() == e.lag.tolist() + + +class TestPredict: + """predict_rating / predict_rating_multi (PlayerRatings predict.rating, + R lines 1056-1133). Pins from the executed oracle + (files/predict_oracle.py); every assert reads wrapper/crate outputs.""" + + def test_elo_branch_pins(self): + import numpy as np + from fast_mlsirm import predict_rating + + p1 = predict_rating([2200.0, 2200.0], [20, 20], [0], [1], gamma=0.0) + assert p1.shape == (1,) + assert p1[0] == 0.5 + p2 = predict_rating([2200.0, 2000.0], [20, 20], [0], [1], gamma=30.0) + assert abs(p2[0] - 0.7898441797581306) < 1e-14 + # per-game gamma array and scalar broadcast agree (crate values) + p2b = predict_rating( + [2200.0, 2000.0], [20, 20], [0, 0], [1, 1], gamma=np.array([30.0, 30.0]) + ) + assert p2b[0] == p2[0] and p2b[1] == p2[0] + + def test_deviation_branch(self): + from fast_mlsirm import predict_rating + + p3 = predict_rating( + [2200.0, 2000.0], [20, 20], [0], [1], + deviations=[50.0, 100.0], gamma=30.0, + ) + assert abs(p3[0] - 0.776912664201114) < 1e-14 + p3z = predict_rating( + [2200.0, 2000.0], [20, 20], [0], [1], + deviations=[0.0, 0.0], gamma=30.0, + ) + p2 = predict_rating([2200.0, 2000.0], [20, 20], [0], [1], gamma=30.0) + assert p3z[0] == p2[0] + + def test_tng_trat_unmatched(self): + import math + from fast_mlsirm import predict_rating + + kept = predict_rating([2200.0, 2000.0], [15, 15], [0], [1], gamma=0.0) + assert abs(kept[0] - 0.7597469266479578) < 1e-14 + dropped = predict_rating([2200.0, 2000.0], [14, 15], [0], [1], gamma=0.0) + assert math.isnan(dropped[0]) + replaced = predict_rating( + [2200.0, 2000.0], [14, 15], [0], [1], gamma=0.0, trat=2000.0 + ) + assert replaced[0] == 0.5 + un = predict_rating([2200.0, 2000.0], [20, 20], [-1], [1], gamma=0.0) + assert math.isnan(un[0]) + un_t = predict_rating( + [2200.0, 2000.0], [20, 20], [-1], [1], gamma=0.0, trat=2000.0 + ) + assert un_t[0] == 0.5 + # stored-NaN + pair trat with deviations + p9c = predict_rating( + [2200.0, 2000.0], [20, 20], [0], [1], + deviations=[float("nan"), 50.0], gamma=30.0, trat=(2200.0, 50.0), + ) + assert abs(p9c[0] - 0.7844611342833985) < 1e-14 + + def test_thresh(self): + import math + from fast_mlsirm import predict_rating + + p = predict_rating( + [2200.0, 2200.0], [20, 20], [0, -1], [1, 1], gamma=0.0, thresh=0.5 + ) + assert p[0] == 1.0 # exact equality maps to 1 (>=) + assert math.isnan(p[1]) + + def test_multi_and_placing(self): + import math + import numpy as np + from fast_mlsirm import predict_rating_multi + + p7 = predict_rating_multi( + [2300.0, 2200.0, 2100.0, 2000.0], [20] * 4, + np.array([[0, 1, 2], [0, 3, -1]]), + ) + assert p7.shape == (2, 3) + assert p7[0].tolist() == [2.5, 0.0, -2.5] + assert p7[1][0] == 3.75 and p7[1][1] == -3.75 and math.isnan(p7[1][2]) + p8 = predict_rating_multi( + [2300.0, 2300.0, 2100.0, 2000.0], [20] * 4, + np.array([[0, 1, 2, -1]]), placing=True, + ) + assert p8[0][0] == 1.0 and p8[0][1] == 1.0 and p8[0][2] == 3.0 + assert math.isnan(p8[0][3]) + + def test_validation(self): + import numpy as np + import pytest + from fast_mlsirm import predict_rating, predict_rating_multi + + with pytest.raises(ValueError, match="2..=10000"): + predict_rating([2200.0], [20], [0], [0], gamma=0.0) + with pytest.raises(ValueError, match="complex"): + predict_rating([2200.0 + 1j, 2000.0], [20, 20], [0], [1]) + with pytest.raises(ValueError, match="complex|bool"): + predict_rating(np.array([True, False]), [20, 20], [0], [1]) + with pytest.raises(ValueError, match="non-numeric"): + predict_rating([2200.0, "x"], [20, 20], [0], [1]) + with pytest.raises(ValueError, match="masked"): + predict_rating( + np.ma.masked_array([2200.0, 2000.0], mask=[False, True]), + [20, 20], [0], [1], + ) + with pytest.raises(ValueError, match="infinit"): + predict_rating([np.inf, 2000.0], [20, 20], [0], [1]) + with pytest.raises(ValueError, match="NaN"): + predict_rating([2200.0, 2000.0], [np.nan, 20], [0], [1]) + with pytest.raises(ValueError, match=">= -1"): + predict_rating([2200.0, 2000.0], [20, 20], [-2], [1]) + with pytest.raises(ValueError, match="self-play"): + predict_rating([2200.0, 2000.0], [20, 20], [1], [1]) + with pytest.raises(ValueError, match="pair"): + predict_rating( + [2200.0, 2000.0], [20, 20], [0], [1], + deviations=[50.0, 50.0], trat=2000.0, + ) + with pytest.raises(ValueError, match="bool"): + predict_rating([2200.0, 2000.0], [20, 20], [0], [1], thresh=True) + with pytest.raises(ValueError, match="2-D"): + predict_rating_multi([2200.0, 2000.0], [20, 20], [0, 1]) + with pytest.raises(ValueError, match="bool"): + predict_rating_multi( + [2200.0, 2000.0], [20, 20], np.array([[0, 1]]), trat=True + ) + # NaN ratings are allowed (R NA), not an error + import math + out = predict_rating([float("nan"), 2000.0], [20, 20], [0], [1], gamma=0.0) + assert math.isnan(out[0]) diff --git a/tests/unit/scaling_tests.rs b/tests/unit/scaling_tests.rs index 1371cc86d..3b847639b 100644 --- a/tests/unit/scaling_tests.rs +++ b/tests/unit/scaling_tests.rs @@ -5615,3 +5615,666 @@ fn fd_mc_500_elo_reduction() { assert_eq!(f.lag, e.lag, "rep {rep}"); } } + +// --------------------------------------------------------------------------- +// predict_rating_two / predict_rating_multi (PlayerRatings predict.rating, +// R lines 1056-1133 READ). Pins from the executed oracle +// (files/predict_oracle.py / predict_oracle_output.txt). +// --------------------------------------------------------------------------- + +/// P1+P2: Elo branch exact half + gamma-sign pin. +/// Asserts read `predict_rating_two` return values. Killing mutants: +/// gamma sign flip (MU1) changes the P2 pin; equal-rating case pins 1/2. +#[test] +fn pr_anchor_elo_gamma_sign() { + let p1 = predict_rating_two( + &[2200.0, 2200.0], + None, + &[20, 20], + &[0], + &[1], + &[0.0], + 15, + None, + None, + ) + .unwrap(); + assert_eq!(p1, vec![0.5]); + let p2 = predict_rating_two( + &[2200.0, 2000.0], + None, + &[20, 20], + &[0], + &[1], + &[30.0], + 15, + None, + None, + ) + .unwrap(); + assert!((p2[0] - 0.7898441797581306).abs() < 1e-14, "{}", p2[0]); + // MU1 (gamma sign flip) would give ~0.7276; assert distance from it. + let mu1 = 1.0 / (1.0 + 10f64.powf((2000.0 - 2200.0 + 30.0) / 400.0)); + assert!((p2[0] - mu1).abs() > 1e-3); +} + +/// P3: deviation-branch shrink pin + dev=0 reduction to the Elo branch. +/// Asserts read `predict_rating_two` values from BOTH branches. Killing +/// mutants: qip3 factor 3->2 (MU2) changes the nonzero-dev pin; the dev=0 +/// crate-vs-crate identity kills stray vec offsets. +#[test] +fn pr_deviation_shrink() { + let p3 = predict_rating_two( + &[2200.0, 2000.0], + Some(&[50.0, 100.0]), + &[20, 20], + &[0], + &[1], + &[30.0], + 15, + None, + None, + ) + .unwrap(); + assert!((p3[0] - 0.776912664201114).abs() < 1e-14, "{}", p3[0]); + let p3z = predict_rating_two( + &[2200.0, 2000.0], + Some(&[0.0, 0.0]), + &[20, 20], + &[0], + &[1], + &[30.0], + 15, + None, + None, + ) + .unwrap(); + let p2 = predict_rating_two( + &[2200.0, 2000.0], + None, + &[20, 20], + &[0], + &[1], + &[30.0], + 15, + None, + None, + ) + .unwrap(); + assert_eq!(p3z[0], p2[0], "dev=0 must reduce to the Elo branch"); +} + +/// P4: tng boundary — games == tng is KEPT (strict <); games < tng is +/// replaced by trat or NaN. Asserts read `predict_rating_two` values. +/// Killing mutant: `<` -> `<=` (MU3) turns the kept pin into NaN. +#[test] +fn pr_tng_boundary() { + let kept = predict_rating_two( + &[2200.0, 2000.0], + None, + &[15, 15], + &[0], + &[1], + &[0.0], + 15, + None, + None, + ) + .unwrap(); + assert!((kept[0] - 0.7597469266479578).abs() < 1e-14, "{}", kept[0]); + let dropped = predict_rating_two( + &[2200.0, 2000.0], + None, + &[14, 15], + &[0], + &[1], + &[0.0], + 15, + None, + None, + ) + .unwrap(); + assert!(dropped[0].is_nan()); + let replaced = predict_rating_two( + &[2200.0, 2000.0], + None, + &[14, 15], + &[0], + &[1], + &[0.0], + 15, + Some((2000.0, 0.0)), + None, + ) + .unwrap(); + assert_eq!(replaced[0], 0.5); +} + +/// P5+P9: unmatched (-1) and matched-but-stored-NaN players — trat +/// replaces ALL extracted NAs; without trat they propagate. Asserts read +/// `predict_rating_two` values (incl. a crate-vs-crate deviation anchor). +#[test] +fn pr_unmatched_and_stored_na_trat() { + let una = predict_rating_two( + &[2200.0, 2000.0], + None, + &[20, 20], + &[-1], + &[1], + &[0.0], + 15, + None, + None, + ) + .unwrap(); + assert!(una[0].is_nan()); + let unb = predict_rating_two( + &[2200.0, 2000.0], + None, + &[20, 20], + &[-1], + &[1], + &[0.0], + 15, + Some((2000.0, 0.0)), + None, + ) + .unwrap(); + assert_eq!(unb[0], 0.5); + let p9a = predict_rating_two( + &[f64::NAN, 2000.0], + None, + &[20, 20], + &[0], + &[1], + &[0.0], + 15, + Some((2000.0, 0.0)), + None, + ) + .unwrap(); + assert_eq!(p9a[0], 0.5); + let p9b = predict_rating_two( + &[f64::NAN, 2000.0], + None, + &[20, 20], + &[0], + &[1], + &[0.0], + 15, + None, + None, + ) + .unwrap(); + assert!(p9b[0].is_nan()); + // stored-NaN deviation replaced by trat.1 == crate value with real dev. + let p9c = predict_rating_two( + &[2200.0, 2000.0], + Some(&[f64::NAN, 50.0]), + &[20, 20], + &[0], + &[1], + &[30.0], + 15, + Some((2200.0, 50.0)), + None, + ) + .unwrap(); + let p9c_ref = predict_rating_two( + &[2200.0, 2000.0], + Some(&[50.0, 50.0]), + &[20, 20], + &[0], + &[1], + &[30.0], + 15, + None, + None, + ) + .unwrap(); + assert_eq!(p9c[0], p9c_ref[0]); + assert!((p9c[0] - 0.7844611342833985).abs() < 1e-14, "{}", p9c[0]); +} + +/// P6: thresh uses >= (exact equality -> 1) and NaN preds stay NaN. +/// Asserts read `predict_rating_two` values. Killing mutant: `>=` -> `>` +/// (MU6) turns the exact-equality 1.0 into 0.0. +#[test] +fn pr_thresh_ge_and_nan() { + let p = predict_rating_two( + &[2200.0, 2200.0], + None, + &[20, 20], + &[0, -1], + &[1, 1], + &[0.0, 0.0], + 15, + None, + Some(0.5), + ) + .unwrap(); + assert_eq!(p[0], 1.0, "pred exactly == thresh must map to 1 (>=)"); + assert!(p[1].is_nan()); + let below = predict_rating_two( + &[2000.0, 2200.0], + None, + &[20, 20], + &[0], + &[1], + &[0.0], + 15, + None, + Some(0.5), + ) + .unwrap(); + assert_eq!(below[0], 0.0); +} + +/// P7: EloM (rat - rowmean)/40 pins, na.rm rowmean, all-NaN row. +/// Asserts read `predict_rating_multi` values. Killing mutant: divisor +/// 40 -> 400 (MU4) changes the exact 2.5/3.75 pins. +#[test] +fn pr_elom_rowmean() { + let p7 = predict_rating_multi( + &[2300.0, 2200.0, 2100.0, 2000.0], + &[20, 20, 20, 20], + &[0, 1, 2, 0, 3, -1], + 2, + 3, + 15, + None, + false, + ) + .unwrap(); + assert_eq!(&p7[0..3], &[2.5, 0.0, -2.5]); + assert_eq!(p7[3], 3.75); + assert_eq!(p7[4], -3.75); + assert!(p7[5].is_nan()); + let allna = predict_rating_multi( + &[2300.0, 2200.0], + &[20, 20], + &[-1, -1], + 1, + 2, + 15, + None, + false, + ) + .unwrap(); + assert!(allna.iter().all(|v| v.is_nan())); + // P9d: stored-NaN seat replaced by trat. + let p9d = predict_rating_multi( + &[f64::NAN, 2200.0, 2100.0], + &[20, 20, 20], + &[0, 1, 2], + 1, + 3, + 15, + Some(2300.0), + false, + ) + .unwrap(); + // tng boundary in the multi branch: games == tng KEPT (strict <). + let boundary = predict_rating_multi( + &[2300.0, 2200.0, 2100.0], + &[15, 20, 20], + &[0, 1, 2], + 1, + 3, + 15, + None, + false, + ) + .unwrap(); + assert_eq!(&boundary[0..3], &[2.5, 0.0, -2.5]); + let dropped = predict_rating_multi( + &[2300.0, 2200.0, 2100.0], + &[14, 20, 20], + &[0, 1, 2], + 1, + 3, + 15, + None, + false, + ) + .unwrap(); + assert!(dropped[0].is_nan()); + assert_eq!(dropped[1], 1.25); + assert_eq!(dropped[2], -1.25); +} + +/// P8: placing ranks — ties share the MINIMUM rank, NaN kept. +/// Asserts read `predict_rating_multi` values. Killing mutant: min -> +/// average/max tie handling (MU5) changes the (1,1,3) pattern. +#[test] +fn pr_placing_min_ties() { + let p8 = predict_rating_multi( + &[2300.0, 2300.0, 2100.0, 2000.0], + &[20, 20, 20, 20], + &[0, 1, 2, -1], + 1, + 4, + 15, + None, + true, + ) + .unwrap(); + assert_eq!(p8[0], 1.0); + assert_eq!(p8[1], 1.0); + assert_eq!(p8[2], 3.0); + assert!(p8[3].is_nan()); +} + +/// Error contract for both functions. Asserts read Err values. +#[test] +fn pr_error_contract() { + let r = |x: Result, String>| x.unwrap_err(); + // player count bounds + assert!(r(predict_rating_two( + &[1500.0], + None, + &[0], + &[0], + &[0], + &[0.0], + 15, + None, + None + )) + .contains("2..=10000")); + // games length mismatch + assert!(r(predict_rating_two( + &[1500.0, 1500.0], + None, + &[0], + &[0], + &[1], + &[0.0], + 15, + None, + None + )) + .contains("games length")); + // empty game rows + assert!(r(predict_rating_two( + &[1500.0, 1500.0], + None, + &[0, 0], + &[], + &[], + &[], + 15, + None, + None + )) + .contains("at least one game")); + // infinite rating rejected (NaN allowed elsewhere) + assert!(r(predict_rating_two( + &[f64::INFINITY, 1500.0], + None, + &[0, 0], + &[0], + &[1], + &[0.0], + 15, + None, + None + )) + .contains("infinite")); + // non-finite gamma + assert!(r(predict_rating_two( + &[1500.0, 1500.0], + None, + &[0, 0], + &[0], + &[1], + &[f64::NAN], + 15, + None, + None + )) + .contains("gamma")); + // non-finite trat + assert!(r(predict_rating_two( + &[1500.0, 1500.0], + None, + &[0, 0], + &[0], + &[1], + &[0.0], + 15, + Some((f64::NAN, 0.0)), + None + )) + .contains("trat")); + // trat deviation checked only when deviations supplied + assert!(predict_rating_two( + &[1500.0, 1500.0], + None, + &[0, 0], + &[0], + &[1], + &[0.0], + 15, + Some((1500.0, f64::NAN)), + None + ) + .is_ok()); + assert!(r(predict_rating_two( + &[1500.0, 1500.0], + Some(&[50.0, 50.0]), + &[0, 0], + &[0], + &[1], + &[0.0], + 15, + Some((1500.0, f64::NAN)), + None + )) + .contains("trat deviation")); + // non-finite thresh + assert!(r(predict_rating_two( + &[1500.0, 1500.0], + None, + &[0, 0], + &[0], + &[1], + &[0.0], + 15, + None, + Some(f64::NAN) + )) + .contains("thresh")); + // out-of-range index (only -1 sentinel allowed) + assert!(r(predict_rating_two( + &[1500.0, 1500.0], + None, + &[0, 0], + &[-2], + &[1], + &[0.0], + 15, + None, + None + )) + .contains("out of range")); + assert!(r(predict_rating_two( + &[1500.0, 1500.0], + None, + &[0, 0], + &[0], + &[2], + &[0.0], + 15, + None, + None + )) + .contains("out of range")); + // self-play + assert!(r(predict_rating_two( + &[1500.0, 1500.0], + None, + &[0, 0], + &[1], + &[1], + &[0.0], + 15, + None, + None + )) + .contains("self-play")); + // deviations length mismatch + assert!(r(predict_rating_two( + &[1500.0, 1500.0], + Some(&[50.0]), + &[0, 0], + &[0], + &[1], + &[0.0], + 15, + None, + None + )) + .contains("deviations length")); + // multi: seat bounds, players length, nr=0, index range, trat + assert!(r(predict_rating_multi( + &[1500.0, 1500.0], + &[0, 0], + &[0], + 1, + 1, + 15, + None, + false + )) + .contains("2..=1000")); + assert!(r(predict_rating_multi( + &[1500.0, 1500.0], + &[0, 0], + &[0, 1, 0], + 1, + 2, + 15, + None, + false + )) + .contains("players length")); + assert!(r(predict_rating_multi( + &[1500.0, 1500.0], + &[0, 0], + &[], + 0, + 2, + 15, + None, + false + )) + .contains("at least one event")); + assert!(r(predict_rating_multi( + &[1500.0, 1500.0], + &[0, 0], + &[0, 5], + 1, + 2, + 15, + None, + false + )) + .contains("out of range")); + assert!(r(predict_rating_multi( + &[1500.0, 1500.0], + &[0, 0], + &[0, 1], + 1, + 2, + 15, + Some(f64::INFINITY), + false + )) + .contains("trat")); +} + +/// MC-500: structural invariants over random inputs. Every assert reads +/// crate return values: complement symmetry pred(w,b)+pred(b,w) ~= 1 at +/// gamma=0 (both branches), preds in (0,1), thresh output consistent with +/// the crate's own unthresholded preds, placing ranks are a valid +/// min-tie ranking of the crate's own EloM preds. +#[test] +#[ignore] +fn pr_mc_500_invariants() { + let mut rng = Lcg(0x5eed_cafe_1234_0001); + for rep in 0..500 { + let n = 3 + (rng.next_f64() * 8.0) as usize; + let ratings: Vec = (0..n).map(|_| 1200.0 + 1600.0 * rng.next_f64()).collect(); + let devs: Vec = (0..n).map(|_| 30.0 + 300.0 * rng.next_f64()).collect(); + let games: Vec = (0..n) + .map(|_| 15 + (rng.next_f64() * 40.0) as u64) + .collect(); + let w = (rng.next_f64() * n as f64) as i64; + let mut b = (rng.next_f64() * n as f64) as i64; + if b == w { + b = (b + 1) % n as i64; + } + for dev_opt in [None, Some(devs.as_slice())] { + let fwd = predict_rating_two( + &ratings, + dev_opt, + &games, + &[w], + &[b], + &[0.0], + 15, + None, + None, + ) + .unwrap(); + let bwd = predict_rating_two( + &ratings, + dev_opt, + &games, + &[b], + &[w], + &[0.0], + 15, + None, + None, + ) + .unwrap(); + assert!(fwd[0] > 0.0 && fwd[0] < 1.0, "rep {rep}"); + assert!( + (fwd[0] + bwd[0] - 1.0).abs() < 1e-12, + "rep {rep}: complement symmetry" + ); + let th = predict_rating_two( + &ratings, + dev_opt, + &games, + &[w], + &[b], + &[0.0], + 15, + None, + Some(0.5), + ) + .unwrap(); + let expect = if fwd[0] >= 0.5 { 1.0 } else { 0.0 }; + assert_eq!(th[0], expect, "rep {rep}: thresh vs crate pred"); + } + // EloM: placing must be the min-tie ranking of the crate's preds. + let np = 3.min(n); + let seats: Vec = (0..np as i64).collect(); + let preds = predict_rating_multi(&ratings, &games, &seats, 1, np, 15, None, false).unwrap(); + let ranks = predict_rating_multi(&ratings, &games, &seats, 1, np, 15, None, true).unwrap(); + for s in 0..np { + let expected = 1.0 + + preds + .iter() + .filter(|u| !u.is_nan() && **u > preds[s]) + .count() as f64; + assert_eq!(ranks[s], expected, "rep {rep} seat {s}"); + } + // Row mean of preds is 0 (crate values; na.rm mean identity). + let sum: f64 = preds.iter().sum(); + assert!(sum.abs() < 1e-9, "rep {rep}: pred row sum {sum}"); + } +} From 4cbba186a6f0aae7629a31abca79bccfda9180e6 Mon Sep 17 00:00:00 2001 From: Seongho Bae Date: Mon, 27 Jul 2026 08:20:04 +0900 Subject: [PATCH 2/2] Fix predict_rating u64 fidelity and multi nr*np overflow MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Impl-review round-1 findings: (1) games/tng went through a float round-trip in the Python wrappers, losing exact integer counts at or above 2^53 and silently shifting the strict games < tng cutoff — new _predict_games_u64/_predict_tng_u64 keep integer inputs lossless and reject float inputs at/above the source dtype's exact-integer bound; (2) predict_rating_multi computed nr*np unchecked, so a wrapping product could pass the length check and panic on indexing — now checked_mul with a checked error. Regression tests at both the Rust error contract and pytest levels. Co-authored-by: Copilot App <223556219+Copilot@users.noreply.github.com> --- crates/mlsirm-core/src/scaling.rs | 9 ++- python/fast_mlsirm/scaling.py | 113 +++++++++++++++++++++++++----- tests/test_paper_features.py | 47 +++++++++++++ tests/unit/scaling_tests.rs | 14 ++++ 4 files changed, 162 insertions(+), 21 deletions(-) diff --git a/crates/mlsirm-core/src/scaling.rs b/crates/mlsirm-core/src/scaling.rs index 1d918dc65..338b7ac2c 100644 --- a/crates/mlsirm-core/src/scaling.rs +++ b/crates/mlsirm-core/src/scaling.rs @@ -3921,11 +3921,14 @@ pub fn predict_rating_multi( "predict_rating_multi: seats per event must be in 2..=1000, got {np}" )); } - if players.len() != nr * np { + let total = nr.checked_mul(np).ok_or_else(|| { + format!("predict_rating_multi: nr * np overflows usize (nr = {nr}, np = {np})") + })?; + if players.len() != total { return Err(format!( "predict_rating_multi: players length {} != nr*np = {}", players.len(), - nr * np + total )); } if let Some(t) = trat { @@ -3940,7 +3943,7 @@ pub fn predict_rating_multi( )); } } - let mut out = vec![f64::NAN; nr * np]; + let mut out = vec![f64::NAN; total]; let mut row_rat = vec![f64::NAN; np]; for i in 0..nr { for s in 0..np { diff --git a/python/fast_mlsirm/scaling.py b/python/fast_mlsirm/scaling.py index 3089faf93..22a97a8f9 100644 --- a/python/fast_mlsirm/scaling.py +++ b/python/fast_mlsirm/scaling.py @@ -1837,6 +1837,87 @@ def _predict_scalar(x, name, fname): return v +def _predict_games_u64(x, fname): + """Validate games counts losslessly into u64. + + Integer-dtype ndarrays (and Python-int sequences, which numpy keeps in + an integer dtype) pass through without a float round-trip, preserving + counts at or above 2**53. Float inputs are bounded by the source + dtype's exact-integer limit (np.finfo(...).nmant + 1 bits: float64 + 2**53, float32 2**24) so a rounded count can never silently shift the + strict games < tng cutoff.""" + name = "games" + if isinstance(x, np.ma.MaskedArray): + raise ValueError(f"{fname}: masked arrays are not supported for {name}") + if isinstance(x, np.ndarray): + raw = x + else: + probe = np.asarray(x, dtype=object) + if np.iscomplexobj(probe): + raise ValueError(f"{fname}: {name} must be real numeric, not complex/bool") + for v in np.ravel(probe): + if v is None or isinstance( + v, (bool, np.bool_, str, bytes, np.datetime64, np.timedelta64) + ): + raise ValueError(f"{fname}: {name} contains a non-numeric value") + try: + raw = np.asarray(x) + except (TypeError, ValueError, OverflowError) as exc: + raise ValueError(f"{fname}: {name} is not numeric: {exc}") from None + if np.iscomplexobj(raw) or raw.dtype.kind == "b": + raise ValueError(f"{fname}: {name} must be real numeric, not complex/bool") + if raw.ndim != 1: + raise ValueError(f"{fname}: {name} must be one-dimensional") + if raw.dtype.kind in "iu": + if raw.dtype.kind == "i" and np.any(raw < 0): + raise ValueError(f"{fname}: {name} must be nonnegative integers") + return raw.astype(np.uint64) + arr = _predict_float_array(x, name, fname, allow_nan=False) + if np.any(arr < 0) or not np.all(arr == np.floor(arr)): + raise ValueError(f"{fname}: {name} must be nonnegative integers") + if raw.dtype.kind == "f": + fidelity = 2.0 ** (np.finfo(raw.dtype).nmant + 1) + else: + fidelity = 2.0**53 + if np.any(arr >= fidelity): + raise ValueError( + f"{fname}: {name} values at or above {int(fidelity)} are not " + "reliably representable in this float dtype; pass games as an " + "integer array" + ) + return arr.astype(np.uint64) + + +def _predict_tng_u64(tng, fname): + """Validate tng losslessly into u64 (no float round-trip for ints).""" + import math + + if isinstance(tng, (bool, np.bool_)): + raise ValueError(f"{fname}: tng must be real numeric, not bool") + if isinstance(tng, (int, np.integer)): + t = int(tng) + if t < 0: + raise ValueError(f"{fname}: tng must be a nonnegative integer") + if t > 2**64 - 1: + raise ValueError(f"{fname}: tng must fit in an unsigned 64-bit integer") + return t + if ( + isinstance(tng, np.ndarray) + and tng.ndim == 0 + and tng.dtype.kind in "iu" + ): + return _predict_tng_u64(int(tng), fname) + v = _predict_scalar(tng, "tng", fname) + if v < 0 or v != math.floor(v): + raise ValueError(f"{fname}: tng must be a nonnegative integer") + if v >= 2.0**53: + raise ValueError( + f"{fname}: tng values at or above 2**53 are not reliably " + "representable as float; pass tng as an int" + ) + return int(v) + + def predict_rating( ratings, games, @@ -1862,8 +1943,6 @@ def predict_rating( predictions to 1 when ``pred >= thresh`` else 0 (NaN kept). ``gamma`` is the per-game (or scalar, broadcast) first-player advantage. """ - import math - from .fitstats import _core_module fname = "predict_rating" @@ -1871,10 +1950,12 @@ def predict_rating( n = rat.shape[0] if n < 2 or n > 10000: raise ValueError(f"{fname}: number of players must be in 2..=10000, got {n}") - g_arr = _predict_float_array(games, "games", fname, allow_nan=False) - if np.any(g_arr < 0) or not np.all(g_arr == np.floor(g_arr)): - raise ValueError(f"{fname}: games must be nonnegative integers") - games_u64 = g_arr.astype(np.uint64) + games_u64 = _predict_games_u64(games, fname) + if games_u64.shape[0] != n: + raise ValueError( + f"{fname}: games must have one entry per player ({n}), " + f"got {games_u64.shape[0]}" + ) w = _predict_int_index_array(white, "white", fname) b = _predict_int_index_array(black, "black", fname) dev = None @@ -1892,9 +1973,7 @@ def predict_rating( raise ValueError( f"{fname}: gamma must be a scalar or length-{ng} array, got {gam.shape}" ) - tng_v = _predict_scalar(tng, "tng", fname) - if tng_v < 0 or tng_v != math.floor(tng_v): - raise ValueError(f"{fname}: tng must be a nonnegative integer") + tng_v = _predict_tng_u64(tng, fname) trat_rating = None trat_deviation = None if trat is not None: @@ -1951,8 +2030,6 @@ def predict_rating_multi( by min-tie ranks of the predictions (rank 1 = highest; NaN kept), matching R's ``rank(-preds, na.last="keep", ties.method="min")``. """ - import math - from .fitstats import _core_module fname = "predict_rating_multi" @@ -1960,10 +2037,12 @@ def predict_rating_multi( n = rat.shape[0] if n < 2 or n > 10000: raise ValueError(f"{fname}: number of players must be in 2..=10000, got {n}") - g_arr = _predict_float_array(games, "games", fname, allow_nan=False) - if np.any(g_arr < 0) or not np.all(g_arr == np.floor(g_arr)): - raise ValueError(f"{fname}: games must be nonnegative integers") - games_u64 = g_arr.astype(np.uint64) + games_u64 = _predict_games_u64(games, fname) + if games_u64.shape[0] != n: + raise ValueError( + f"{fname}: games must have one entry per player ({n}), " + f"got {games_u64.shape[0]}" + ) if isinstance(players, np.ma.MaskedArray): raise ValueError(f"{fname}: masked arrays are not supported for players") p_raw = players if isinstance(players, np.ndarray) else np.asarray(players, dtype=object) @@ -1975,9 +2054,7 @@ def predict_rating_multi( raise ValueError( f"{fname}: seats per event must be in 2..=1000, got {np_seats}" ) - tng_v = _predict_scalar(tng, "tng", fname) - if tng_v < 0 or tng_v != math.floor(tng_v): - raise ValueError(f"{fname}: tng must be a nonnegative integer") + tng_v = _predict_tng_u64(tng, fname) trat_v = None if trat is None else _predict_scalar(trat, "trat", fname) if not isinstance(placing, (bool, np.bool_)): raise ValueError(f"{fname}: placing must be a bool") diff --git a/tests/test_paper_features.py b/tests/test_paper_features.py index 8ea26dfae..47706adf0 100644 --- a/tests/test_paper_features.py +++ b/tests/test_paper_features.py @@ -10201,3 +10201,50 @@ def test_validation(self): import math out = predict_rating([float("nan"), 2000.0], [20, 20], [0], [1], gamma=0.0) assert math.isnan(out[0]) + + def test_u64_fidelity_and_overflow(self): + import math + import numpy as np + import pytest + from fast_mlsirm import predict_rating, predict_rating_multi + + big = 2**53 + # int-dtype games keep exact counts above 2**53: games == 2**53 is + # strictly below tng == 2**53 + 1, so both players are unrated -> NaN. + out = predict_rating( + [2200.0, 2000.0], [big, big], [0], [1], gamma=0.0, tng=big + 1 + ) + assert math.isnan(out[0]) + outm = predict_rating_multi( + [2300.0, 2200.0], [big, big], np.array([[0, 1]]), tng=big + 1 + ) + assert math.isnan(outm[0][0]) and math.isnan(outm[0][1]) + # float games at/above the dtype exact-integer bound are rejected + with pytest.raises(ValueError, match='integer array'): + predict_rating( + [2200.0, 2000.0], np.array([float(big), 20.0]), [0], [1] + ) + with pytest.raises(ValueError, match='integer array'): + predict_rating_multi( + [2200.0, 2000.0], np.array([float(big), 20.0]), np.array([[0, 1]]) + ) + # float tng at/above 2**53 rejected; huge int tng still exact + with pytest.raises(ValueError, match='pass tng as an int'): + predict_rating([2200.0, 2000.0], [20, 20], [0], [1], tng=float(big)) + kept = predict_rating( + [2200.0, 2000.0], [big + 2, big + 2], [0], [1], gamma=0.0, tng=big + 1 + ) + assert abs(kept[0] - 0.7597469266479578) < 1e-14 + # core-level nr*np overflow is a checked ValueError, not a panic + import fast_mlsirm._core as core + with pytest.raises(ValueError, match='overflows'): + core.predict_rating_multi( + np.array([1.0, 2.0]), + np.array([1, 1], dtype=np.uint64), + np.array([], dtype=np.int64), + 1 << 61, + 8, + 15, + None, + False, + ) diff --git a/tests/unit/scaling_tests.rs b/tests/unit/scaling_tests.rs index 3b847639b..5cd041dbd 100644 --- a/tests/unit/scaling_tests.rs +++ b/tests/unit/scaling_tests.rs @@ -6159,6 +6159,20 @@ fn pr_error_contract() { false )) .contains("players length")); + // nr*np usize overflow must be a checked error, not a wrap + panic + // (kills: replacing checked_mul with wrapping `nr * np`, which lets an + // empty players slice pass the length check and then index OOB). + assert!(r(predict_rating_multi( + &[1500.0, 1500.0], + &[0, 0], + &[], + 1usize << 61, + 8, + 15, + None, + false + )) + .contains("overflows")); assert!(r(predict_rating_multi( &[1500.0, 1500.0], &[0, 0],