ferrolearn-preprocess 0.5.0

//! Divergence audit for `ferrolearn-preprocess/src/feature_selection.rs`
//! (`VarianceThreshold`, `SelectKBest`, `anova_f_scores`) against
//! scikit-learn 1.5.2:
//! - `VarianceThreshold` — `sklearn/feature_selection/_variance_threshold.py`
//! - `SelectKBest` / `_get_support_mask` / `_clean_nans` —
//!   `sklearn/feature_selection/_univariate_selection.py`
//!
//! Tracking: #1424.
//!
//! Oracle values (sklearn 1.5.2, LIVE) — hard-coded per R-CHAR-3. Generated by:
//! ```text
//! python3 -c "import numpy as np
//! from sklearn.feature_selection import SelectKBest, f_classif, VarianceThreshold
//! import warnings; warnings.simplefilter('ignore')
//! X=np.array([[1.,1.,5.],[2.,2.,5.],[7.,7.,5.],[8.,8.,5.]]); y=np.array([0,0,1,1])
//! sk=SelectKBest(f_classif,k=1).fit(X,y)
//! print(sk.scores_.tolist())                                  # [72.0, 72.0, nan]
//! print(np.flatnonzero(sk.get_support()).tolist())            # k=1 -> [1]
//! print(np.flatnonzero(SelectKBest(f_classif,k=2).fit(X,y).get_support()).tolist())  # k=2 -> [0,1]
//! print(np.flatnonzero(SelectKBest(f_classif,k=10).fit(X,y).get_support()).tolist()) # k=10 -> [0,1,2]
//! Xd=np.array([[1.,1.,1.],[2.,3.,5.],[8.,7.,2.],[9.,8.,6.]])
//! print(SelectKBest(f_classif,k=2).fit(Xd,y).scores_.tolist()) # [98.0, 24.2, 0.125]
//! print(np.flatnonzero(SelectKBest(f_classif,k=2).fit(Xd,y).get_support()).tolist()) # [0,1]
//! Xv=np.array([[1.,5.,2.],[2.,5.,2.],[3.,5.,2.]])
//! vt=VarianceThreshold(0.0).fit(Xv)
//! print(vt.variances_.tolist())                               # [0.6667,0,0]
//! print(np.flatnonzero(vt.get_support()).tolist())            # [0]"
//! ```

use ferrolearn_core::traits::Fit;
use ferrolearn_preprocess::feature_selection::{ScoreFunc, SelectKBest, VarianceThreshold};
use ndarray::{Array1, Array2, array};

// ---------------------------------------------------------------------------
// Oracle constants (sklearn 1.5.2, LIVE — hard-coded per R-CHAR-3)
// ---------------------------------------------------------------------------

/// DIV-A/B fixture: cols 0,1 identical (tie at F=72.0), col 2 constant (F=NaN).
/// y = [0,0,1,1].
const TIE_SCORES_INFORMATIVE: [f64; 2] = [72.0, 72.0]; // sklearn scores_[0], scores_[1]
/// k=1 support: sklearn keeps the HIGHER index of the tie (argsort mergesort[-k:]).
const TIE_K1_SUPPORT: [usize; 1] = [1];
/// k=2 support: the two informative tied cols; the constant col 2 (NaN -> finfo.min)
/// ranks LAST and is excluded.
const TIE_K2_SUPPORT: [usize; 2] = [0, 1];
/// k=10 (> n_features=3) support: sklearn warns and keeps ALL features.
const K_OVER_NFEAT_SUPPORT: [usize; 3] = [0, 1, 2];

/// Distinct-finite-score fixture (no ties, no constants): F = [98.0, 24.2, 0.125].
const DISTINCT_SCORES: [f64; 3] = [98.0, 24.2, 0.125];
/// k=2 on the distinct fixture: top two by score.
const DISTINCT_K2_SUPPORT: [usize; 2] = [0, 1];

/// VarianceThreshold population variances for Xv = [[1,5,2],[2,5,2],[3,5,2]].
const VT_VARIANCES: [f64; 3] = [0.666_666_666_666_666_6, 0.0, 0.0];

// ===========================================================================
// DIV-A — SelectKBest k-boundary tie-break: higher vs lower index
// ===========================================================================
//
// sklearn `_get_support_mask` (`_univariate_selection.py:794`):
//   `mask[np.argsort(scores, kind="mergesort")[-self.k:]] = 1`
// argsort is ascending + stable; taking the last k means that on a TIE at the
// k-boundary the element appearing LATER in the array (higher index) wins.
// ferrolearn (`feature_selection.rs:489-495`) sorts DESCENDING with a
// `.then(a.cmp(&b))` tie-break that keeps the LOWER index.
//
// Fixture: cols 0 and 1 are identical (F=72.0 tie); col 2 is constant (NaN).
// k=1: sklearn keeps col [1]; ferrolearn keeps col [0].

/// Divergence: `SelectKBest::fit` diverges from
/// `sklearn/feature_selection/_univariate_selection.py:794`
/// (`mask[np.argsort(scores, kind="mergesort")[-self.k:]] = 1`) on a
/// k-boundary tie. sklearn keeps the HIGHER index ([1]); ferrolearn keeps the
/// LOWER index ([0]).
/// Tracking: #1424.
#[test]
fn divergence_selectkbest_tiebreak_keeps_higher_index() {
    let x: Array2<f64> = array![[1., 1., 5.], [2., 2., 5.], [7., 7., 5.], [8., 8., 5.]];
    let y: Array1<usize> = array![0, 0, 1, 1];
    let sel = SelectKBest::<f64>::new(1, ScoreFunc::FClassif);
    let fitted = match sel.fit(&x, &y) {
        Ok(f) => f,
        Err(e) => panic!("fit failed: {e:?}"),
    };
    // sklearn k=1 support is [1] (higher index of the F=72.0 tie).
    assert_eq!(
        fitted.selected_indices(),
        &TIE_K1_SUPPORT[..],
        "sklearn argsort-mergesort[-1:] keeps the HIGHER index of the tie"
    );
}

// ===========================================================================
// DIV-B — SelectKBest constant-feature score: NaN (ranks last) vs +inf (ranks first)
// ===========================================================================
//
// In sklearn a constant feature yields F = 0/0 = NaN; `_clean_nans`
// (`_univariate_selection.py:31-32`) maps NaN -> `np.finfo(dtype).min`, so it
// ranks LAST and is never selected unless k >= n_features.
// ferrolearn `anova_f_scores` (`feature_selection.rs:424-426`) returns
// `F::infinity()` whenever `ms_within == 0`, including the constant case
// (ms_between == 0 too), ranking the constant feature FIRST -> wrongly selected.
//
// Fixture: col 2 is constant. k=2: sklearn selects the two informative cols
// [0,1]; ferrolearn selects the constant col 2 (+inf) plus one informative col.

/// Divergence: `anova_f_scores` (`feature_selection.rs:424-426`) returns +inf
/// for a CONSTANT feature, whereas sklearn yields NaN ->
/// `_clean_nans` (`_univariate_selection.py:31-32`) -> `finfo.min`, ranking it
/// LAST. With k=2 sklearn selects the two informative cols [0,1]; ferrolearn
/// wrongly includes the constant col 2.
/// Tracking: #1424.
#[test]
fn divergence_selectkbest_constant_feature_ranks_last() {
    let x: Array2<f64> = array![[1., 1., 5.], [2., 2., 5.], [7., 7., 5.], [8., 8., 5.]];
    let y: Array1<usize> = array![0, 0, 1, 1];
    let sel = SelectKBest::<f64>::new(2, ScoreFunc::FClassif);
    let fitted = match sel.fit(&x, &y) {
        Ok(f) => f,
        Err(e) => panic!("fit failed: {e:?}"),
    };
    // sklearn k=2 support is [0,1] — the constant col 2 (NaN) ranks last.
    assert_eq!(
        fitted.selected_indices(),
        &TIE_K2_SUPPORT[..],
        "constant col 2 (NaN->finfo.min in sklearn) must NOT be selected at k=2"
    );
}

// ===========================================================================
// DIV-C — SelectKBest k > n_features: sklearn warns + keeps all; ferrolearn errors
// ===========================================================================
//
// sklearn `SelectKBest._check_params` (`_univariate_selection.py:774-779`)
// WARNS (does NOT raise) when k > n_features and returns all features.
// LIVE ORACLE: SelectKBest(f_classif, k=10).fit(X_3feat, y) succeeds with
// support [0,1,2]. ferrolearn `fit` (`feature_selection.rs:470-478`) returns
// `Err(InvalidParameter)`. This is a divergence.
// Tracking: separate blocker (see report) + #1424.

/// Divergence: `SelectKBest::fit` (`feature_selection.rs:470-478`) returns
/// `Err(InvalidParameter)` when `k > n_features`, whereas sklearn
/// `_check_params` (`_univariate_selection.py:774-779`) only WARNS and keeps
/// all features. LIVE oracle: k=10 on a 3-feature X -> support [0,1,2].
/// Tracking: #1424.
#[test]
fn divergence_selectkbest_k_over_nfeatures_keeps_all() {
    let x: Array2<f64> = array![[1., 1., 5.], [2., 2., 5.], [7., 7., 5.], [8., 8., 5.]];
    let y: Array1<usize> = array![0, 0, 1, 1];
    let sel = SelectKBest::<f64>::new(10, ScoreFunc::FClassif);
    let fitted = match sel.fit(&x, &y) {
        Ok(f) => f,
        Err(e) => panic!("sklearn keeps all features when k>n_features; ferrolearn errored: {e:?}"),
    };
    assert_eq!(
        fitted.selected_indices(),
        &K_OVER_NFEAT_SUPPORT[..],
        "sklearn warns and keeps ALL features when k>n_features"
    );
}

// ===========================================================================
// GREEN-GUARDS — behavior that MUST match sklearn (un-ignored)
// ===========================================================================

/// Green-guard: ANOVA F-scores on finite non-constant features match sklearn
/// `scores_` for the tied informative cols 0,1 (F=72.0 each).
#[test]
fn guard_selectkbest_scores_match_sklearn_finite() {
    let x: Array2<f64> = array![[1., 1., 5.], [2., 2., 5.], [7., 7., 5.], [8., 8., 5.]];
    let y: Array1<usize> = array![0, 0, 1, 1];
    let sel = SelectKBest::<f64>::new(2, ScoreFunc::FClassif);
    let fitted = match sel.fit(&x, &y) {
        Ok(f) => f,
        Err(e) => panic!("fit failed: {e:?}"),
    };
    let scores = fitted.scores();
    assert!(
        (scores[0] - TIE_SCORES_INFORMATIVE[0]).abs() < 1e-9,
        "col0 F={} expected {}",
        scores[0],
        TIE_SCORES_INFORMATIVE[0]
    );
    assert!(
        (scores[1] - TIE_SCORES_INFORMATIVE[1]).abs() < 1e-9,
        "col1 F={} expected {}",
        scores[1],
        TIE_SCORES_INFORMATIVE[1]
    );
}

/// Green-guard: normal top-k (no ties, no constants, distinct F-scores).
/// sklearn k=2 support is [0,1]; scores [98.0, 24.2, 0.125].
#[test]
fn guard_selectkbest_distinct_topk_matches() {
    let x: Array2<f64> = array![[1., 1., 1.], [2., 3., 5.], [8., 7., 2.], [9., 8., 6.]];
    let y: Array1<usize> = array![0, 0, 1, 1];
    let sel = SelectKBest::<f64>::new(2, ScoreFunc::FClassif);
    let fitted = match sel.fit(&x, &y) {
        Ok(f) => f,
        Err(e) => panic!("fit failed: {e:?}"),
    };
    let scores = fitted.scores();
    for j in 0..3 {
        assert!(
            (scores[j] - DISTINCT_SCORES[j]).abs() < 1e-9,
            "col{j} F={} expected {}",
            scores[j],
            DISTINCT_SCORES[j]
        );
    }
    assert_eq!(fitted.selected_indices(), &DISTINCT_K2_SUPPORT[..]);
}

/// Green-guard: a perfect-separator single feature (ss_between>0, ms_within=0)
/// yields F=+inf in BOTH sklearn and ferrolearn — the +inf is CORRECT here
/// (NOT a divergence). sklearn score is inf; ferrolearn keeps it at k=1.
#[test]
fn guard_selectkbest_perfect_separator_is_inf_and_kept() {
    let x: Array2<f64> = array![[1.], [1.], [9.], [9.]];
    let y: Array1<usize> = array![0, 0, 1, 1];
    let sel = SelectKBest::<f64>::new(1, ScoreFunc::FClassif);
    let fitted = match sel.fit(&x, &y) {
        Ok(f) => f,
        Err(e) => panic!("fit failed: {e:?}"),
    };
    // sklearn f_classif returns inf for this perfect separator.
    assert!(
        fitted.scores()[0].is_infinite(),
        "perfect separator F must be +inf (matches sklearn inf)"
    );
    assert_eq!(fitted.selected_indices(), &[0usize][..]);
}

/// Green-guard: VarianceThreshold(0.0) — population variances + support match
/// sklearn. Xv col1,col2 constant -> dropped; col0 kept.
#[test]
fn guard_variance_threshold_zero_drops_constants() {
    let x: Array2<f64> = array![[1., 5., 2.], [2., 5., 2.], [3., 5., 2.]];
    let sel = VarianceThreshold::<f64>::new(0.0);
    let fitted = match sel.fit(&x, &()) {
        Ok(f) => f,
        Err(e) => panic!("fit failed: {e:?}"),
    };
    let v = fitted.variances();
    for j in 0..3 {
        assert!(
            (v[j] - VT_VARIANCES[j]).abs() < 1e-9,
            "var col{j}={} expected {}",
            v[j],
            VT_VARIANCES[j]
        );
    }
    assert_eq!(fitted.selected_indices(), &[0usize][..]);
}

/// Green-guard: VarianceThreshold custom threshold 0.5.
/// Xc col0 var=0.6667 (>0.5 kept), col1 var=66.67 (kept) -> support [0,1].
#[test]
fn guard_variance_threshold_custom_threshold() {
    let x: Array2<f64> = array![[1., 10.], [2., 20.], [3., 30.]];
    let sel = VarianceThreshold::<f64>::new(0.5);
    let fitted = match sel.fit(&x, &()) {
        Ok(f) => f,
        Err(e) => panic!("fit failed: {e:?}"),
    };
    // col0 var = 0.6667 > 0.5; col1 var = 66.67 > 0.5 — both kept.
    assert_eq!(fitted.selected_indices(), &[0usize, 1][..]);
    let v = fitted.variances();
    assert!((v[0] - 0.666_666_666_666_666_6).abs() < 1e-9);
    assert!((v[1] - 66.666_666_666_666_67).abs() < 1e-9);
}

/// Green-guard: VarianceThreshold rejects a negative threshold.
#[test]
fn guard_variance_threshold_negative_threshold_errors() {
    let x: Array2<f64> = array![[1.], [2.]];
    let sel = VarianceThreshold::<f64>::new(-0.1);
    assert!(sel.fit(&x, &()).is_err());
}

/// Green-guard: VarianceThreshold rejects zero rows.
#[test]
fn guard_variance_threshold_zero_rows_errors() {
    let x: Array2<f64> = Array2::zeros((0, 2));
    let sel = VarianceThreshold::<f64>::new(0.0);
    assert!(sel.fit(&x, &()).is_err());
}

/// Green-guard: SelectKBest rejects zero rows.
#[test]
fn guard_selectkbest_zero_rows_errors() {
    let x: Array2<f64> = Array2::zeros((0, 3));
    let y: Array1<usize> = Array1::zeros(0);
    let sel = SelectKBest::<f64>::new(1, ScoreFunc::FClassif);
    assert!(sel.fit(&x, &y).is_err());
}

/// Green-guard: SelectKBest rejects a y-length mismatch.
#[test]
fn guard_selectkbest_y_length_mismatch_errors() {
    let x: Array2<f64> = array![[1., 2.], [3., 4.]];
    let y: Array1<usize> = array![0]; // wrong length
    let sel = SelectKBest::<f64>::new(1, ScoreFunc::FClassif);
    assert!(sel.fit(&x, &y).is_err());
}

/// Green-guard: f32 path — VarianceThreshold(0.0) drops the constant col.
#[test]
fn guard_variance_threshold_f32() {
    let x: Array2<f32> = array![[1.0f32, 5.0], [2.0, 5.0], [3.0, 5.0]];
    let sel = VarianceThreshold::<f32>::new(0.0f32);
    let fitted = match sel.fit(&x, &()) {
        Ok(f) => f,
        Err(e) => panic!("fit failed: {e:?}"),
    };
    assert_eq!(fitted.selected_indices(), &[0usize][..]);
}

// ===========================================================================
// RE-AUDIT (#1424) — post-fix faithfulness probes vs LIVE sklearn 1.5.2.
//
// Oracle values hard-coded per R-CHAR-3, generated LIVE (warnings suppressed):
//   python3 -c "import numpy as np; from sklearn.feature_selection import \
//   SelectKBest, f_classif; import warnings; warnings.simplefilter('ignore'); \
//   ...; print(np.flatnonzero(SelectKBest(f_classif,k=K).fit(X,y) \
//   .get_support()).tolist())"
// `_clean_nans` (`_univariate_selection.py:24-33`) maps NaN -> finfo.min;
// `_get_support_mask` (`_univariate_selection.py:794`) does
// `mask[argsort(scores, kind="mergesort")[-k:]] = 1` (ascending + STABLE; a
// k-boundary tie keeps the HIGHER index).
//
// The `assert!(false, ...)` in fit-failure arms is a deliberate test-failure
// signal; each such fn carries a per-fn allow for the constant-assertion lint.
// ===========================================================================

/// RE-AUDIT (a): four IDENTICAL features (all F=72.0), k=2.
/// LIVE oracle: support == [2, 3] (the two HIGHEST indices of the 4-way tie).
#[test]
#[allow(
    clippy::assertions_on_constants,
    reason = "assert!(false,..) is the test-failure signal in the fit-Err arm"
)]
fn reaudit_a_four_way_tie_k2_keeps_two_highest() {
    let x: Array2<f64> = array![
        [1., 1., 1., 1.],
        [2., 2., 2., 2.],
        [7., 7., 7., 7.],
        [8., 8., 8., 8.]
    ];
    let y: Array1<usize> = array![0, 0, 1, 1];
    let sel = SelectKBest::<f64>::new(2, ScoreFunc::FClassif);
    let fitted = match sel.fit(&x, &y) {
        Ok(f) => f,
        Err(e) => {
            assert!(false, "fit failed: {e:?}");
            return;
        }
    };
    assert_eq!(
        fitted.selected_indices(),
        &[2usize, 3][..],
        "4-way tie k=2 must keep the two HIGHEST indices (sklearn argsort-mergesort[-2:])"
    );
}

/// RE-AUDIT (b): 5 features — cols 0,1 CONSTANT (NaN score), cols 2,3,4
/// informative with distinct F = [98.0, 24.2, 0.125].
/// LIVE oracle: k=3 -> [2, 3, 4] (the 3 informative; constants rank last).
/// k=4 -> [1, 2, 3, 4] (forced to include ONE constant: argsort of cleaned
/// scores is [0,1,4,3,2]; last-4 keeps the HIGHER-index constant col 1).
#[test]
#[allow(
    clippy::assertions_on_constants,
    reason = "assert!(false,..) is the test-failure signal in the fit-Err arm"
)]
fn reaudit_b_two_constants_plus_informative_k3_k4() {
    let x: Array2<f64> = array![
        [9., 9., 1., 1., 1.],
        [9., 9., 2., 3., 5.],
        [9., 9., 8., 7., 2.],
        [9., 9., 9., 8., 6.]
    ];
    let y: Array1<usize> = array![0, 0, 1, 1];

    let sel3 = SelectKBest::<f64>::new(3, ScoreFunc::FClassif);
    let f3 = match sel3.fit(&x, &y) {
        Ok(f) => f,
        Err(e) => {
            assert!(false, "k=3 fit failed: {e:?}");
            return;
        }
    };
    assert_eq!(
        f3.selected_indices(),
        &[2usize, 3, 4][..],
        "k=3 must pick the 3 informative cols; both constants (NaN->finfo.min) rank last"
    );

    let sel4 = SelectKBest::<f64>::new(4, ScoreFunc::FClassif);
    let f4 = match sel4.fit(&x, &y) {
        Ok(f) => f,
        Err(e) => {
            assert!(false, "k=4 fit failed: {e:?}");
            return;
        }
    };
    assert_eq!(
        f4.selected_indices(),
        &[1usize, 2, 3, 4][..],
        "k=4 forces one constant; sklearn keeps the HIGHER-index constant (col 1, not col 0)"
    );
}

/// RE-AUDIT (c): k-boundary triple — k==n_features (all), k==0 (none),
/// k==1 distinct (single argmax). Fixture F = [98.0, 24.2, 0.125].
/// LIVE oracle: k=3 -> [0,1,2]; k=0 -> []; k=1 -> [0].
#[test]
#[allow(
    clippy::assertions_on_constants,
    reason = "assert!(false,..) is the test-failure signal in the fit-Err arm"
)]
fn reaudit_c_k_boundaries_all_none_argmax() {
    let x: Array2<f64> = array![[1., 1., 1.], [2., 3., 5.], [8., 7., 2.], [9., 8., 6.]];
    let y: Array1<usize> = array![0, 0, 1, 1];

    let f_all = match SelectKBest::<f64>::new(3, ScoreFunc::FClassif).fit(&x, &y) {
        Ok(f) => f,
        Err(e) => {
            assert!(false, "k=all fit failed: {e:?}");
            return;
        }
    };
    assert_eq!(
        f_all.selected_indices(),
        &[0usize, 1, 2][..],
        "k==n_features selects all"
    );

    let f_none = match SelectKBest::<f64>::new(0, ScoreFunc::FClassif).fit(&x, &y) {
        Ok(f) => f,
        Err(e) => {
            assert!(false, "k=0 fit failed: {e:?}");
            return;
        }
    };
    assert!(
        f_none.selected_indices().is_empty(),
        "k==0 selects NONE (sklearn _get_support_mask k==0 -> all-false), got {:?}",
        f_none.selected_indices()
    );

    let f_one = match SelectKBest::<f64>::new(1, ScoreFunc::FClassif).fit(&x, &y) {
        Ok(f) => f,
        Err(e) => {
            assert!(false, "k=1 fit failed: {e:?}");
            return;
        }
    };
    assert_eq!(
        f_one.selected_indices(),
        &[0usize][..],
        "k==1 picks the single argmax (col 0, F=98)"
    );
}

/// RE-AUDIT (d): k=7 on 4 features (k > n_features).
/// LIVE oracle: sklearn warns and keeps ALL: support == [0,1,2,3].
#[test]
#[allow(
    clippy::assertions_on_constants,
    reason = "assert!(false,..) is the test-failure signal in the fit-Err arm"
)]
fn reaudit_d_k_over_nfeatures_clamps_keep_all() {
    let x: Array2<f64> = array![
        [1., 1., 1., 1.],
        [2., 3., 5., 2.],
        [8., 7., 2., 9.],
        [9., 8., 6., 1.]
    ];
    let y: Array1<usize> = array![0, 0, 1, 1];
    let sel = SelectKBest::<f64>::new(7, ScoreFunc::FClassif);
    let fitted = match sel.fit(&x, &y) {
        Ok(f) => f,
        Err(e) => {
            assert!(false, "k>n_features must keep all (warn), not error: {e:?}");
            return;
        }
    };
    assert_eq!(
        fitted.selected_indices(),
        &[0usize, 1, 2, 3][..],
        "k=7 on 4 features -> clamp + keep all (sklearn warn+keep-all)"
    );
}

/// RE-AUDIT (e): MIXED — cols 0,1 TIE (F=72.0), col 2 CONSTANT (NaN),
/// cols 3,4 distinct (F=24.2, 0.05). cleaned argsort = [2,4,3,0,1].
/// LIVE oracle: k=3 -> [0,1,3]; k=4 -> [0,1,3,4].
#[test]
#[allow(
    clippy::assertions_on_constants,
    reason = "assert!(false,..) is the test-failure signal in the fit-Err arm"
)]
fn reaudit_e_mixed_tie_constant_distinct() {
    let x: Array2<f64> = array![
        [1., 1., 5., 1., 1.],
        [2., 2., 5., 3., 9.],
        [7., 7., 5., 7., 2.],
        [8., 8., 5., 8., 6.]
    ];
    let y: Array1<usize> = array![0, 0, 1, 1];

    let f3 = match SelectKBest::<f64>::new(3, ScoreFunc::FClassif).fit(&x, &y) {
        Ok(f) => f,
        Err(e) => {
            assert!(false, "k=3 fit failed: {e:?}");
            return;
        }
    };
    assert_eq!(
        f3.selected_indices(),
        &[0usize, 1, 3][..],
        "mixed k=3: the two tied cols (0,1) + the higher distinct (3); constant col 2 excluded"
    );

    let f4 = match SelectKBest::<f64>::new(4, ScoreFunc::FClassif).fit(&x, &y) {
        Ok(f) => f,
        Err(e) => {
            assert!(false, "k=4 fit failed: {e:?}");
            return;
        }
    };
    assert_eq!(
        f4.selected_indices(),
        &[0usize, 1, 3, 4][..],
        "mixed k=4: adds the lower distinct col 4; constant col 2 still last"
    );
}

/// RE-AUDIT (f): real-ish 8x4, 3 classes, no exact ties, k=2.
/// LIVE oracle scores_ = [38.6639..., 0.79534..., 287.7576..., 43.3094...];
/// k=2 support == [2, 3] (the two highest F).
#[test]
#[allow(
    clippy::assertions_on_constants,
    reason = "assert!(false,..) is the test-failure signal in the fit-Err arm"
)]
fn reaudit_f_realish_8x4_3classes_k2() {
    let x: Array2<f64> = array![
        [5.1, 3.5, 1.4, 0.2],
        [4.9, 3.0, 1.4, 0.2],
        [6.2, 2.2, 4.5, 1.5],
        [5.9, 3.2, 4.8, 1.8],
        [7.3, 2.9, 6.3, 1.8],
        [6.7, 3.3, 5.7, 2.5],
        [5.0, 3.4, 1.5, 0.2],
        [6.4, 3.2, 4.5, 1.5]
    ];
    let y: Array1<usize> = array![0, 0, 1, 1, 2, 2, 0, 1];
    let sel = SelectKBest::<f64>::new(2, ScoreFunc::FClassif);
    let fitted = match sel.fit(&x, &y) {
        Ok(f) => f,
        Err(e) => {
            assert!(false, "fit failed: {e:?}");
            return;
        }
    };
    let scores = fitted.scores();
    let expected = [
        38.663_903_061_219_145,
        0.795_347_744_360_936_9,
        287.757_601_351_380_97,
        43.309_426_229_507_73,
    ];
    for (j, &want) in expected.iter().enumerate() {
        assert!(
            (scores[j] - want).abs() < 1e-9,
            "col{j} F={} expected {want}",
            scores[j]
        );
    }
    assert_eq!(
        fitted.selected_indices(),
        &[2usize, 3][..],
        "realish k=2 -> the two highest-F cols [2,3]"
    );
}

/// RE-AUDIT (g): scores() parity with a constant feature. sklearn scores_ is
/// NaN at the constant col (col 2) and finite elsewhere (cols 0,1 tie F=72.0,
/// cols 3,4 distinct). Verify ferrolearn scores() is NaN at col 2 and matches
/// finite scores to 1e-9.
#[test]
#[allow(
    clippy::assertions_on_constants,
    reason = "assert!(false,..) is the test-failure signal in the fit-Err arm"
)]
fn reaudit_g_scores_nan_at_constant_finite_elsewhere() {
    let x: Array2<f64> = array![
        [1., 1., 5., 1., 1.],
        [2., 2., 5., 3., 9.],
        [7., 7., 5., 7., 2.],
        [8., 8., 5., 8., 6.]
    ];
    let y: Array1<usize> = array![0, 0, 1, 1];
    // LIVE oracle scores_ = [72.0, 72.0, nan, 24.2, 0.05].
    let sel = SelectKBest::<f64>::new(3, ScoreFunc::FClassif);
    let fitted = match sel.fit(&x, &y) {
        Ok(f) => f,
        Err(e) => {
            assert!(false, "fit failed: {e:?}");
            return;
        }
    };
    let s = fitted.scores();
    assert!(
        s[2].is_nan(),
        "constant col 2 score must be NaN (matches sklearn scores_[2]), got {}",
        s[2]
    );
    let finite = [(0usize, 72.0), (1, 72.0), (3, 24.2), (4, 0.05)];
    for (j, want) in finite {
        assert!(
            (s[j] - want).abs() < 1e-9,
            "col{j} F={} expected {want}",
            s[j]
        );
    }
}

/// RE-AUDIT (h): f32 path with a 4-way tie, k=2.
/// LIVE oracle (f32): support == [2, 3] (the two highest indices of the tie).
#[test]
#[allow(
    clippy::assertions_on_constants,
    reason = "assert!(false,..) is the test-failure signal in the fit-Err arm"
)]
fn reaudit_h_f32_tie_k2_keeps_two_highest() {
    let x: Array2<f32> = array![
        [1.0f32, 1.0, 1.0, 1.0],
        [2.0, 2.0, 2.0, 2.0],
        [7.0, 7.0, 7.0, 7.0],
        [8.0, 8.0, 8.0, 8.0]
    ];
    let y: Array1<usize> = array![0, 0, 1, 1];
    let sel = SelectKBest::<f32>::new(2, ScoreFunc::FClassif);
    let fitted = match sel.fit(&x, &y) {
        Ok(f) => f,
        Err(e) => {
            assert!(false, "f32 fit failed: {e:?}");
            return;
        }
    };
    assert_eq!(
        fitted.selected_indices(),
        &[2usize, 3][..],
        "f32 4-way tie k=2 must keep the two HIGHEST indices"
    );
}