model-selection-rs 0.1.0

//! The whole point of nested CV, as an explicit test (plan Milestone 6 DoD):
//! a *naive* "tune and report on the same folds" procedure is optimistically
//! biased, while nested CV is not.
//!
//! Construction: labels and every candidate model's predictions are pure noise,
//! uncorrelated with each other. With many candidates, whichever one happens to
//! fit the evaluation folds best will score above chance *on those folds* — that
//! is the optimism. Nested CV selects on inner folds and reports on untouched
//! outer folds, so it recovers the true ~50% chance level.

use ndarray::{Array1, Array2};

use model_selection_rs::evaluate::{cross_validate, nested_cross_validate, BoxedScorer};
use model_selection_rs::scoring::Accuracy;
use model_selection_rs::splitters::KFold;

/// Deterministic pseudo-random bit from a (seed, index) pair (splitmix64-ish).
fn bit(seed: u64, i: u64) -> f64 {
    let mut z = seed
        .wrapping_mul(0x9E37_79B9_7F4A_7C15)
        .wrapping_add(i.wrapping_mul(0xBF58_476D_1CE4_E5B9));
    z ^= z >> 27;
    z = z.wrapping_mul(0x94D0_49BB_1331_11EB);
    z ^= z >> 31;
    (z & 1) as f64
}

const N: usize = 160;
const N_CANDIDATES: u64 = 80;
const Y_SEED: u64 = 12345;

/// Candidate model `p`: predictions are noise seeded by `p`, read off the row-id
/// column so they are a genuine function of the (subset-selected) feature matrix.
fn candidate_fit(
    p: &u64,
    _x: &Array2<f64>,
    _y: &Array1<f64>,
) -> impl Fn(&Array2<f64>) -> Array1<f64> {
    let p = *p;
    move |xq: &Array2<f64>| xq.column(0).mapv(|id| bit(p, id as u64))
}

fn make_data() -> (Array2<f64>, Array1<f64>) {
    // x[i,0] is the row id, carried through index selection so predictions stay
    // tied to sample identity.
    let x = Array2::from_shape_fn((N, 1), |(i, _)| i as f64);
    let y = Array1::from_shape_fn(N, |i| bit(Y_SEED, i as u64));
    (x, y)
}

/// Tune: pick the candidate with the best inner-CV accuracy on the given data.
fn tune_best_candidate(x: &Array2<f64>, y: &Array1<f64>, inner: &KFold) -> u64 {
    let mut best = (f64::NEG_INFINITY, 0u64);
    for p in 0..N_CANDIDATES {
        let scorers: Vec<BoxedScorer> = vec![Box::new(Accuracy)];
        let res = cross_validate(
            inner,
            x,
            y,
            move |xt, yt| candidate_fit(&p, xt, yt),
            &scorers,
            false,
        )
        .unwrap();
        let acc = res.mean_test_score(0);
        if acc > best.0 {
            best = (acc, p);
        }
    }
    best.1
}

#[test]
fn naive_selection_is_optimistic_nested_is_not() {
    let (x, y) = make_data();

    // ---- Naive estimate: tune over candidates on the full data, then report
    // that same best cross-validated score. This double-uses every fold.
    let cv = KFold::new(5).unwrap().with_shuffle(7);
    let mut naive_best = f64::NEG_INFINITY;
    for p in 0..N_CANDIDATES {
        let scorers: Vec<BoxedScorer> = vec![Box::new(Accuracy)];
        let res = cross_validate(
            &cv,
            &x,
            &y,
            move |xt, yt| candidate_fit(&p, xt, yt),
            &scorers,
            false,
        )
        .unwrap();
        naive_best = naive_best.max(res.mean_test_score(0));
    }

    // ---- Nested estimate: tune on inner folds, score on untouched outer folds.
    let outer = KFold::new(5).unwrap().with_shuffle(7);
    let inner = KFold::new(4).unwrap().with_shuffle(13);
    let nested = nested_cross_validate(
        &outer,
        &inner,
        &x,
        &y,
        tune_best_candidate,
        candidate_fit,
        &Accuracy,
    )
    .unwrap();
    let nested_mean = nested.mean_score();

    // The labels are noise: the honest estimate should sit near chance...
    assert!(
        (0.40..=0.60).contains(&nested_mean),
        "nested estimate {nested_mean} should be near chance (0.5)"
    );
    // ...while the naive best-of-many is inflated well above chance...
    assert!(
        naive_best > 0.55,
        "naive best {naive_best} should be optimistically above chance"
    );
    // ...and strictly more optimistic than the nested estimate.
    assert!(
        naive_best > nested_mean + 0.02,
        "naive {naive_best} should exceed nested {nested_mean} — the optimism bias"
    );

    // The per-outer-fold selections are also surfaced as a diagnostic.
    assert_eq!(nested.selected_params.len(), 5);
}