1use std::collections::{BTreeMap, HashMap};
77
78use serde::{Deserialize, Serialize};
79
80use super::summary::{Apply, Check as Logged, Gate, Label, Summary};
81use super::{ARMS, Route, Unit};
82use crate::conscience::{self, Challenge};
83use crate::memory::Rng;
84
85pub const ALPHA: f64 = 0.10;
92pub const DELTA: f64 = 0.05;
95pub const DRAWS: usize = conscience::P_DRAWS;
98pub const SOLVE_MAX: usize = 16;
103
104pub const PRIOR_SD: f64 = 2.5;
108pub const ASK_PRIOR: (f64, f64) = (0.0, 0.35);
111pub const NOT_PRIOR: (f64, f64) = (-0.3, 0.25);
116pub const DILUTION_PRIOR: (f64, f64) = (0.0, 0.1);
120pub const ARD_START: f64 = 1.0;
123pub const LEVEL_START: f64 = 0.5;
125pub const SD_FLOOR: f64 = 0.01;
128pub const HALF_LIVES: [Option<u32>; 5] = [None, Some(64), Some(32), Some(16), Some(8)];
131pub const FIT_SESSIONS: usize = 256;
134pub const MAX_IPW: f64 = 10.0;
137pub const PSI_FLOOR: f64 = 1e-3;
141pub const POOL_STEPS: [i32; 5] = [-2, -1, 0, 1, 2];
144pub const LATENCY_MIN: usize = 20;
147pub const LATENCY_NEAR: f64 = 0.9;
154
155pub const ROUTE_LABELS_MIN: usize = 5;
164
165const TOL: f64 = 1e-8;
167const MAX_NEWTON: usize = 50;
169const MAX_HALVINGS: usize = 30;
172const MAX_EB: usize = 50;
174const EB_TOL: f64 = 1e-4;
176const PREQ_SESSIONS: usize = 16;
179const CS_SESSIONS: f64 = 30.0;
182const CS_MIN: usize = 10;
185const PRICE_STEPS: usize = 60;
187
188const B0: usize = 0;
190const G2: usize = 1;
191const G3: usize = 2;
192const RANK: usize = 3;
193const STOP: usize = 4;
194const THREE: usize = 5;
195const ASK: usize = 6;
196const NOT: usize = 7;
197const ASK_M: usize = 8;
198const FREE_M: usize = 9;
199const ASK_STOP: usize = 10;
200const ASK_G3: usize = 11;
201const NOT_G: usize = 12;
202const POOL: usize = 13;
203const THREE_STOP: usize = 14;
204pub const P: usize = 15;
207pub const NAMES: [&str; P] = [
209 "intercept",
210 "grade2",
211 "grade3",
212 "rank",
213 "stop",
214 "three",
215 "ask",
216 "not",
217 "ask_dilution",
218 "free_dilution",
219 "ask_stop",
220 "ask_grade3",
221 "not_grade",
222 "pool",
223 "three_stop",
224];
225const SESSION: usize = 0;
227const LESSON: usize = 1;
228const JUDGE: usize = 2;
229const FAMILIES: usize = 3;
230
231#[derive(Debug, Clone, Copy)]
233struct Coef {
234 mean: f64,
235 sd: f64,
236 ard: bool,
237 nonpos: bool,
238}
239
240fn coefs() -> [Coef; P] {
243 let fixed = |(mean, sd): (f64, f64), nonpos| Coef { mean, sd, ard: false, nonpos };
244 let ard = |nonpos| Coef { mean: 0.0, sd: ARD_START, ard: true, nonpos };
245 let mut c = [fixed((0.0, PRIOR_SD), false); P];
246 c[RANK] = ard(false);
247 c[ASK] = fixed(ASK_PRIOR, false);
248 c[NOT] = fixed(NOT_PRIOR, false);
249 c[ASK_M] = fixed(DILUTION_PRIOR, true);
250 c[FREE_M] = fixed(DILUTION_PRIOR, true);
251 for i in [ASK_STOP, ASK_G3, NOT_G] {
252 c[i] = ard(false);
253 }
254 c[POOL] = ard(true);
255 c
256}
257
258const PRUNED_KEEPS: [usize; 7] = [B0, G2, G3, ASK, NOT, ASK_M, FREE_M];
261
262fn pruned() -> [Coef; P] {
265 let mut c = coefs();
266 for (k, x) in c.iter_mut().enumerate() {
267 if !PRUNED_KEEPS.contains(&k) {
268 *x = Coef { mean: 0.0, sd: SD_FLOOR, ard: false, nonpos: false };
269 }
270 }
271 c
272}
273
274#[derive(Debug, Clone, Copy, PartialEq)]
276pub struct Feat {
277 pub grade: u8,
278 pub rank: usize,
279 pub stop: bool,
281 pub three: bool,
283 pub pool: f64,
285}
286
287impl Feat {
288 pub fn of(units: &[Unit], stop: bool, tokens: u32) -> Vec<Feat> {
290 let three = units.iter().any(|u| u.grade == 3);
291 let pool = pool_x(tokens);
292 units.iter().map(|u| Feat { grade: u.grade, rank: u.rank, stop, three, pool }).collect()
293 }
294}
295
296fn pool_x(tokens: u32) -> f64 {
298 if tokens == 0 { 0.0 } else { (f64::from(tokens) / crate::whole::TOKENS as f64).ln() }
299}
300
301fn design(f: &Feat, route: Route, m: usize) -> [f64; P] {
303 let mut x = [0.0; P];
304 let g = f.grade;
305 x[B0] = 1.0;
306 x[G2] = f64::from(u8::from(g == 2));
307 x[G3] = f64::from(u8::from(g == 3));
308 x[RANK] = (1.0 + f.rank as f64).ln();
309 x[STOP] = f64::from(u8::from(f.stop));
310 x[THREE] = f64::from(u8::from(f.three));
311 x[THREE_STOP] = f64::from(u8::from(f.stop && f.three));
312 x[POOL] = f.pool;
313 let ln_m = (m.max(1) as f64).ln();
314 match route {
315 Route::Ask => {
316 x[ASK] = 1.0;
317 x[ASK_M] = ln_m;
318 x[ASK_STOP] = x[STOP];
319 x[ASK_G3] = x[G3];
320 }
321 Route::Free => x[FREE_M] = ln_m,
322 Route::Placebo | Route::Not => {
323 x[NOT] = 1.0;
324 x[NOT_G] = f64::from(g) - 2.0;
325 }
326 }
327 x
328}
329
330fn dot<const K: usize>(a: &[f64; K], b: &[f64; K]) -> f64 {
331 a.iter().zip(b).map(|(x, y)| x * y).sum()
332}
333
334fn sigmoid(x: f64) -> f64 {
335 if x >= 0.0 {
336 1.0 / (1.0 + (-x).exp())
337 } else {
338 let e = x.exp();
339 e / (1.0 + e)
340 }
341}
342
343fn log1pexp(x: f64) -> f64 {
345 if x > 0.0 { x + (-x).exp().ln_1p() } else { x.exp().ln_1p() }
346}
347
348fn cholesky(a: &mut [f64], n: usize) -> bool {
351 for j in 0..n {
352 let d = a[j * n + j] - (0..j).map(|k| a[j * n + k] * a[j * n + k]).sum::<f64>();
353 if d.is_nan() || d <= 0.0 || !d.is_finite() {
354 return false;
355 }
356 let d = d.sqrt();
357 a[j * n + j] = d;
358 for i in j + 1..n {
359 let s = a[i * n + j] - (0..j).map(|k| a[i * n + k] * a[j * n + k]).sum::<f64>();
360 a[i * n + j] = s / d;
361 }
362 for i in 0..j {
363 a[i * n + j] = 0.0;
364 }
365 }
366 true
367}
368
369fn chol_solve(l: &[f64], n: usize, b: &mut [f64]) {
371 for i in 0..n {
372 let s = b[i] - (0..i).map(|k| l[i * n + k] * b[k]).sum::<f64>();
373 b[i] = s / l[i * n + i];
374 }
375 for i in (0..n).rev() {
376 let s = b[i] - (i + 1..n).map(|k| l[k * n + i] * b[k]).sum::<f64>();
377 b[i] = s / l[i * n + i];
378 }
379}
380
381fn chol_inverse(l: &[f64], n: usize) -> Vec<f64> {
383 let mut inv = vec![0.0; n * n];
384 let mut e = vec![0.0; n];
385 for j in 0..n {
386 e.iter_mut().for_each(|x| *x = 0.0);
387 e[j] = 1.0;
388 chol_solve(l, n, &mut e);
389 for (i, x) in e.iter().enumerate() {
390 inv[i * n + j] = *x;
391 }
392 }
393 inv
394}
395
396fn normal(rng: &mut Rng) -> f64 {
398 let (u, v) = (rng.unit(), rng.unit());
399 (-2.0 * u.ln()).sqrt() * (std::f64::consts::TAU * v).cos()
400}
401
402fn gamma(shape: f64, rng: &mut Rng) -> f64 {
404 if shape < 1.0 {
405 return gamma(shape + 1.0, rng) * rng.unit().powf(1.0 / shape);
406 }
407 let d = shape - 1.0 / 3.0;
408 let c = 1.0 / (9.0 * d).sqrt();
409 for _ in 0..64 {
412 let z = normal(rng);
413 let v = (1.0 + c * z).powi(3);
414 if v <= 0.0 {
415 continue;
416 }
417 if rng.unit().ln() < 0.5 * z * z + d - d * v + d * v.ln() {
418 return d * v;
419 }
420 }
421 d
422}
423
424fn beta_draw(a: f64, b: f64, rng: &mut Rng) -> f64 {
426 let (x, y) = (gamma(a, rng), gamma(b, rng));
427 if x + y > 0.0 { x / (x + y) } else { a / (a + b) }
428}
429
430fn hash(parts: &[&[u8]]) -> u64 {
432 let mut h = blake3::Hasher::new();
433 h.update(b"borg-select-policy-v1\0");
434 for p in parts {
435 h.update(p);
436 h.update(b"\0");
437 }
438 let mut first = [0u8; 8];
439 first.copy_from_slice(&h.finalize().as_bytes()[..8]);
440 u64::from_le_bytes(first)
441}
442
443fn mix(seed: u64, i: u64) -> u64 {
445 Rng::new(seed ^ i.wrapping_add(1).wrapping_mul(0x9e37_79b9_7f4a_7c15)).next_u64()
446}
447
448#[derive(Debug, Clone)]
450struct Obs<const K: usize, const F: usize> {
451 x: [f64; K],
452 y: f64,
453 w: f64,
454 lv: [Option<usize>; F],
455}
456
457#[derive(Debug, Clone)]
459struct Est<const K: usize, const F: usize> {
460 beta: [f64; K],
461 u: [Vec<f64>; F],
462 sd: [f64; K],
464 spread: [f64; F],
465 cov: [[f64; K]; K],
467 var_u: [Vec<f64>; F],
468 converged: bool,
469}
470
471impl<const K: usize, const F: usize> Est<K, F> {
472 fn eta(&self, o: &Obs<K, F>) -> f64 {
473 dot(&self.beta, &o.x) + (0..F).filter_map(|f| o.lv[f].map(|l| self.u[f][l])).sum::<f64>()
474 }
475}
476
477fn map<const K: usize, const F: usize>(
482 obs: &[Obs<K, F>],
483 levels: [usize; F],
484 spec: &[Coef; K],
485 sd: &[f64; K],
486 spread: &[f64; F],
487 warm: Option<&Est<K, F>>,
488) -> Option<Est<K, F>> {
489 let elim = (0..F).filter(|&f| levels[f] > 0).max_by_key(|&f| (levels[f], std::cmp::Reverse(f)));
490 let mut off = [0usize; F];
491 let mut dim = K;
492 for f in 0..F {
493 if Some(f) != elim {
494 off[f] = dim;
495 dim += levels[f];
496 }
497 }
498 let members: Vec<Vec<usize>> = match elim {
499 Some(e) => {
500 let mut m = vec![Vec::new(); levels[e]];
501 for (i, o) in obs.iter().enumerate() {
502 if let Some(l) = o.lv[e] {
503 m[l].push(i);
504 }
505 }
506 m
507 }
508 None => Vec::new(),
509 };
510 let mut beta: [f64; K] = std::array::from_fn(|k| warm.map_or(spec[k].mean, |w| w.beta[k]));
511 for k in 0..K {
512 if spec[k].nonpos {
513 beta[k] = beta[k].min(0.0);
514 }
515 }
516 let mut u: [Vec<f64>; F] = std::array::from_fn(|f| match warm {
517 Some(w) if w.u[f].len() == levels[f] => w.u[f].clone(),
518 _ => vec![0.0; levels[f]],
519 });
520 let prec: [f64; K] = std::array::from_fn(|k| 1.0 / (sd[k] * sd[k]));
521 let uprec: [f64; F] = std::array::from_fn(|f| 1.0 / (spread[f] * spread[f]));
522 let objective = |beta: &[f64; K], u: &[Vec<f64>; F]| -> f64 {
523 let mut ll = 0.0;
524 for o in obs {
525 let eta = dot(beta, &o.x) + (0..F).filter_map(|f| o.lv[f].map(|l| u[f][l])).sum::<f64>();
526 ll += o.w * (o.y * eta - log1pexp(eta));
527 }
528 let pb: f64 = (0..K).map(|k| (beta[k] - spec[k].mean).powi(2) * prec[k]).sum();
529 let pu: f64 = (0..F).map(|f| u[f].iter().map(|x| x * x).sum::<f64>() * uprec[f]).sum();
530 ll - 0.5 * (pb + pu)
531 };
532 let mut converged = false;
533 let mut s = vec![0.0; dim * dim];
534 let mut d_e: Vec<f64> = Vec::new();
535 let mut b_e: Vec<Vec<(usize, f64)>> = Vec::new();
536 let mut current = objective(&beta, &u);
537 for _ in 0..MAX_NEWTON {
538 let mut g = vec![0.0; dim];
540 s.iter_mut().for_each(|x| *x = 0.0);
541 let mut r = vec![0.0; obs.len()];
542 let mut wt = vec![0.0; obs.len()];
543 let mut dense: Vec<(usize, f64)> = Vec::with_capacity(K + F);
544 for (i, o) in obs.iter().enumerate() {
545 let eta = dot(&beta, &o.x) + (0..F).filter_map(|f| o.lv[f].map(|l| u[f][l])).sum::<f64>();
546 let mu = sigmoid(eta);
547 r[i] = o.w * (o.y - mu);
548 wt[i] = o.w * mu * (1.0 - mu);
549 dense.clear();
550 dense.extend((0..K).filter(|&k| o.x[k] != 0.0).map(|k| (k, o.x[k])));
551 dense.extend((0..F).filter(|&f| Some(f) != elim).filter_map(|f| o.lv[f].map(|l| (off[f] + l, 1.0))));
552 for &(a, xa) in &dense {
553 g[a] += r[i] * xa;
554 for &(b, xb) in &dense {
555 s[a * dim + b] += wt[i] * xa * xb;
556 }
557 }
558 }
559 for k in 0..K {
560 g[k] -= (beta[k] - spec[k].mean) * prec[k];
561 s[k * dim + k] += prec[k];
562 }
563 for f in (0..F).filter(|&f| Some(f) != elim) {
564 for l in 0..levels[f] {
565 g[off[f] + l] -= u[f][l] * uprec[f];
566 s[(off[f] + l) * dim + off[f] + l] += uprec[f];
567 }
568 }
569 let mut g_e: Vec<f64> = Vec::new();
571 d_e.clear();
572 b_e.clear();
573 if let Some(e) = elim {
574 for (l, rows) in members.iter().enumerate() {
575 let mut d = uprec[e];
576 let mut gl = -u[e][l] * uprec[e];
577 let mut b: Vec<(usize, f64)> = Vec::new();
578 for &i in rows {
579 d += wt[i];
580 gl += r[i];
581 let o = &obs[i];
582 for k in 0..K {
583 add_sparse(&mut b, k, wt[i] * o.x[k]);
584 }
585 for f in (0..F).filter(|&f| f != e) {
586 if let Some(lv) = o.lv[f] {
587 add_sparse(&mut b, off[f] + lv, wt[i]);
588 }
589 }
590 }
591 for &(a, va) in &b {
592 for &(c, vc) in &b {
593 s[a * dim + c] -= va * vc / d;
594 }
595 g[a] -= va * gl / d;
596 }
597 d_e.push(d);
598 g_e.push(gl);
599 b_e.push(b);
600 }
601 }
602 let mut held: Vec<usize> = (0..K).filter(|&k| spec[k].nonpos && beta[k] >= 0.0 && g[k] > 0.0).collect();
605 let step = loop {
606 let mut delta = vec![0.0; dim];
607 for &k in &held {
608 delta[k] = -beta[k];
609 }
610 let free: Vec<usize> = (0..dim).filter(|i| !held.contains(i)).collect();
611 let n = free.len();
612 let mut a = vec![0.0; n * n];
613 let mut rhs = vec![0.0; n];
614 for (fi, &i) in free.iter().enumerate() {
615 rhs[fi] = g[i] - held.iter().map(|&h| s[i * dim + h] * delta[h]).sum::<f64>();
616 for (fj, &j) in free.iter().enumerate() {
617 a[fi * n + fj] = s[i * dim + j];
618 }
619 }
620 if !cholesky(&mut a, n) {
621 return None;
622 }
623 chol_solve(&a, n, &mut rhs);
624 for (fi, &i) in free.iter().enumerate() {
625 delta[i] = rhs[fi];
626 }
627 let out: Vec<usize> = (0..K).filter(|&k| spec[k].nonpos && !held.contains(&k) && beta[k] + delta[k] > 0.0).collect();
628 if out.is_empty() {
629 break delta;
630 }
631 held.extend(out);
632 };
633 let delta_e: Vec<f64> = (0..d_e.len()).map(|l| (g_e[l] - b_e[l].iter().map(|&(a, v)| v * step[a]).sum::<f64>()) / d_e[l]).collect();
634 let largest = step.iter().chain(&delta_e).fold(0.0f64, |m, x| m.max(x.abs()));
635 let mut t = 1.0;
637 let mut accepted = false;
638 for _ in 0..=MAX_HALVINGS {
639 let nb: [f64; K] = std::array::from_fn(|k| {
640 let x = beta[k] + t * step[k];
641 if spec[k].nonpos { x.min(0.0) } else { x }
642 });
643 let mut nu = u.clone();
644 for f in 0..F {
645 for l in 0..levels[f] {
646 nu[f][l] += t * if Some(f) == elim { delta_e[l] } else { step[off[f] + l] };
647 }
648 }
649 let next = objective(&nb, &nu);
650 if next >= current - 1e-12 * current.abs().max(1.0) {
652 (beta, u, current) = (nb, nu, next);
653 accepted = true;
654 break;
655 }
656 t *= 0.5;
657 }
658 if !accepted || largest * t < TOL {
659 converged = accepted || largest < TOL;
660 break;
661 }
662 }
663 if !cholesky(&mut s, dim) {
665 return None;
666 }
667 let inv = chol_inverse(&s, dim);
668 let cov: [[f64; K]; K] = std::array::from_fn(|a| std::array::from_fn(|b| inv[a * dim + b]));
669 let var_u: [Vec<f64>; F] = std::array::from_fn(|f| {
670 if Some(f) == elim {
671 (0..levels[f])
672 .map(|l| {
673 let b = &b_e[l];
674 let q: f64 = b.iter().map(|&(a, va)| b.iter().map(|&(c, vc)| va * inv[a * dim + c] * vc).sum::<f64>()).sum();
675 1.0 / d_e[l] + q / (d_e[l] * d_e[l])
676 })
677 .collect()
678 } else {
679 (0..levels[f]).map(|l| inv[(off[f] + l) * dim + off[f] + l]).collect()
680 }
681 });
682 Some(Est { beta, u, sd: *sd, spread: *spread, cov, var_u, converged })
683}
684
685fn add_sparse(b: &mut Vec<(usize, f64)>, i: usize, v: f64) {
687 if v == 0.0 {
688 return;
689 }
690 match b.iter_mut().find(|(j, _)| *j == i) {
691 Some(e) => e.1 += v,
692 None => b.push((i, v)),
693 }
694}
695
696fn eb<const K: usize, const F: usize>(obs: &[Obs<K, F>], levels: [usize; F], spec: &[Coef; K], warm: Option<&Est<K, F>>) -> Option<Est<K, F>> {
700 let mut sd: [f64; K] = std::array::from_fn(|k| warm.filter(|_| spec[k].ard).map_or(spec[k].sd, |w| w.sd[k]));
701 let mut spread: [f64; F] = std::array::from_fn(|f| warm.map_or(LEVEL_START, |w| w.spread[f]));
702 let mut est = map(obs, levels, spec, &sd, &spread, warm)?;
703 for _ in 0..MAX_EB {
704 let nsd: [f64; K] = std::array::from_fn(|k| match spec[k].ard {
705 true => ((est.beta[k] - spec[k].mean).powi(2) + est.cov[k][k]).sqrt().clamp(SD_FLOOR, PRIOR_SD),
706 false => sd[k],
707 });
708 let nspread: [f64; F] = std::array::from_fn(|f| match levels[f] {
709 0 => spread[f],
710 n => (est.u[f].iter().zip(&est.var_u[f]).map(|(x, v)| x * x + v).sum::<f64>() / n as f64).sqrt().clamp(SD_FLOOR, PRIOR_SD),
711 });
712 let moved = sd.iter().zip(&nsd).chain(spread.iter().zip(&nspread)).any(|(a, b)| (a - b).abs() > EB_TOL * a.abs().max(SD_FLOOR));
713 (sd, spread) = (nsd, nspread);
714 est = map(obs, levels, spec, &sd, &spread, Some(&est))?;
715 if !moved {
716 break;
717 }
718 }
719 Some(est)
720}
721
722#[derive(Debug, Clone, Default, PartialEq, Serialize, Deserialize)]
725pub struct Rates {
726 pub free_stop: f64,
728 pub free_live: f64,
730 pub later: f64,
733 pub wait: f64,
735 pub ask_live: f64,
741}
742
743#[derive(Debug, Clone, Default, PartialEq, Serialize, Deserialize)]
745pub struct PoolFit {
746 #[serde(default, skip_serializing_if = "Vec::is_empty")]
749 pub yields: Vec<(f64, f64)>,
750 #[serde(default, skip_serializing_if = "is_default")]
752 pub per_lesson: f64,
753 #[serde(default, skip_serializing_if = "is_default")]
756 pub latency: (f64, f64, f64),
757 #[serde(default, skip_serializing_if = "is_default")]
758 pub latency_n: u32,
759 #[serde(default, skip_serializing_if = "is_default")]
763 pub latency_reach: u32,
764 #[serde(default, skip_serializing_if = "is_default")]
767 pub eta: f64,
768 #[serde(default, skip_serializing_if = "is_default")]
769 pub slope: (f64, f64),
770 #[serde(default, skip_serializing_if = "is_default")]
771 pub not: (f64, f64),
772}
773
774#[derive(Debug, Clone, Default, PartialEq, Serialize, Deserialize)]
776pub struct Model {
777 pub beta: Vec<f64>,
780 pub chol: Vec<f64>,
781 pub sd: Vec<f64>,
783 pub spread: [f64; FAMILIES],
785 pub sessions: BTreeMap<String, (f64, f64)>,
787 pub lessons: BTreeMap<String, (f64, f64)>,
788 pub judges: BTreeMap<String, (f64, f64)>,
789 #[serde(default, skip_serializing_if = "String::is_empty")]
791 pub judge: String,
792 #[serde(default, skip_serializing_if = "std::ops::Not::not")]
793 pub converged: bool,
794}
795
796impl Model {
797 #[must_use]
801 pub fn current(&self) -> bool {
802 self.beta.len() == P && self.chol.len() == P * P && self.sd.len() == P
803 }
804}
805
806#[derive(Debug, Clone, Default, PartialEq, Serialize, Deserialize)]
809pub struct Fit {
810 #[serde(default, skip_serializing_if = "String::is_empty")]
813 pub key: String,
814 #[serde(default, skip_serializing_if = "Option::is_none")]
816 pub model: Option<Model>,
817 #[serde(default, skip_serializing_if = "Option::is_none")]
819 pub half_life: Option<u32>,
820 #[serde(default, skip_serializing_if = "is_default")]
821 pub rates: Rates,
822 #[serde(default, skip_serializing_if = "is_default")]
823 pub pool: PoolFit,
824 #[serde(default, skip_serializing_if = "Option::is_none")]
827 pub budget: Option<f64>,
828 #[serde(default, skip_serializing_if = "is_default")]
829 pub lambda: f64,
830 #[serde(default, skip_serializing_if = "std::ops::Not::not")]
833 pub sensitive: bool,
834 #[serde(default, skip_serializing_if = "is_default")]
836 pub labels: u32,
837 #[serde(default, skip_serializing_if = "is_default")]
838 pub dropped: u32,
839}
840
841fn is_default<T: Default + PartialEq>(v: &T) -> bool {
842 *v == T::default()
843}
844
845pub(super) fn decay(age: u32, h: Option<u32>) -> f64 {
847 h.map_or(1.0, |h| 0.5f64.powf(f64::from(age) / f64::from(h)))
848}
849
850struct Data {
853 obs: Vec<Obs<P, FAMILIES>>,
854 tick: Vec<u32>,
857 base_w: Vec<f64>,
858 unclear: Vec<bool>,
859 levels: [usize; FAMILIES],
860 names: [Vec<String>; FAMILIES],
861 lv_of: HashMap<usize, [Option<usize>; FAMILIES]>,
863 judged_p: HashMap<usize, f64>,
864 dropped: u32,
865}
866
867fn level_of(names: &mut [Vec<String>; FAMILIES], index: &mut [HashMap<String, usize>; FAMILIES], tagged: &(usize, String)) -> usize {
869 let (f, name) = tagged;
870 *index[*f].entry(name.clone()).or_insert_with(|| {
871 names[*f].push(name.clone());
872 names[*f].len() - 1
873 })
874}
875
876fn model_route(a: &Apply) -> Option<Route> {
879 match a.route {
880 Route::Ask | Route::Free if a.delivered != super::summary::Delivered::Dropped => Some(a.route),
881 Route::Placebo => Some(Route::Not),
882 _ => None,
883 }
884}
885
886fn applies_of<'a>(s: &'a Summary, c: &Logged) -> &'a [Apply] {
888 s.applies.get(c.first as usize..c.first as usize + usize::from(c.n)).unwrap_or(&[])
889}
890
891fn logged_had(s: &Summary, c: &Logged) -> Vec<String> {
894 applies_of(s, c).iter().filter(|a| a.nofree).filter_map(|a| s.lessons.get(a.lesson as usize).cloned()).collect()
895}
896
897fn check_view<'a>(s: &'a Summary, c: &Logged) -> (Vec<Unit<'a>>, Vec<Feat>, Vec<Route>) {
899 let applies = applies_of(s, c);
900 let units: Vec<Unit> =
901 applies.iter().map(|a| Unit { lesson: s.lessons.get(a.lesson as usize).map_or("", String::as_str), grade: a.grade, rank: a.rank as usize }).collect();
902 let feats = Feat::of(&units, c.stop, c.tokens);
903 let routes = applies.iter().map(|a| if a.route == Route::Placebo { Route::Not } else { a.route }).collect();
904 (units, feats, routes)
905}
906
907impl Data {
908 fn of(s: &Summary, zero: bool) -> Data {
911 let now = s.sessions.len() as u32;
912 let from = now.saturating_sub(FIT_SESSIONS as u32);
913 let mut names: [Vec<String>; FAMILIES] = Default::default();
914 let mut index: [HashMap<String, usize>; FAMILIES] = Default::default();
915 let judged_p = judged(s);
916 let mut d = Data {
917 obs: Vec::new(),
918 tick: Vec::new(),
919 base_w: Vec::new(),
920 unclear: Vec::new(),
921 levels: [0; FAMILIES],
922 names: Default::default(),
923 lv_of: HashMap::new(),
924 judged_p: HashMap::new(),
925 dropped: 0,
926 };
927 let names_of = |c: &Logged, a: &Apply| -> [Option<(usize, String)>; FAMILIES] {
928 [
929 Some((SESSION, s.sessions.get(c.tick as usize).map_or("", |x| x.id.as_str()).to_string())),
930 s.lessons.get(a.lesson as usize).map(|n| (LESSON, n.clone())),
931 a.judge.and_then(|j| s.judges.get(j as usize).cloned()).map(|n| (JUDGE, n)),
932 ]
933 };
934 for c in &s.checks {
935 let (_, feats, _) = check_view(s, c);
936 let (m_q, m_f) = (usize::from(c.m_q), usize::from(c.m_f));
937 for (j, f) in feats.iter().enumerate() {
938 let i = c.first as usize + j;
939 let a = &s.applies[i];
940 let Some(route) = model_route(a) else { continue };
941 let y = match a.label {
942 Label::Changed => 1.0,
943 Label::Ignored => 0.0,
944 Label::Unclear if zero => 0.0,
945 _ => continue,
946 };
947 if c.tick < from {
948 d.dropped += 1;
949 continue;
950 }
951 let m = match route {
952 Route::Ask => m_q,
953 Route::Free => m_f,
954 _ => 1,
955 };
956 let p = judged_p.get(&i).copied().unwrap_or(1.0);
957 let w = 1.0 / p.max(1.0 / MAX_IPW);
958 let lv = names_of(c, a).map(|n| n.map(|n| level_of(&mut names, &mut index, &n)));
959 d.obs.push(Obs { x: design(f, route, m), y, w, lv });
960 d.tick.push(c.tick);
961 d.base_w.push(w);
962 d.unclear.push(a.label == Label::Unclear);
963 }
964 }
965 for c in &s.checks {
967 for i in c.first as usize..(c.first as usize + c.n as usize).min(s.applies.len()) {
968 let n = names_of(c, &s.applies[i]);
969 d.lv_of.insert(i, std::array::from_fn(|f| n[f].as_ref().and_then(|(_, n)| index[f].get(n).copied())));
970 }
971 }
972 d.levels = std::array::from_fn(|f| names[f].len());
973 d.names = names;
974 d.judged_p = judged_p;
975 d
976 }
977
978 fn without_levels(mut self) -> Data {
980 for o in &mut self.obs {
981 o.lv = [None; FAMILIES];
982 }
983 for lv in self.lv_of.values_mut() {
984 *lv = [None; FAMILIES];
985 }
986 self.levels = [0; FAMILIES];
987 self.names = Default::default();
988 self
989 }
990
991 fn weighted(&self, now: u32, h: Option<u32>) -> Vec<Obs<P, FAMILIES>> {
993 self.obs.iter().zip(&self.tick).zip(&self.base_w).map(|((o, t), w)| Obs { w: w * decay(now.saturating_sub(*t), h), ..o.clone() }).collect()
994 }
995}
996
997pub(crate) fn judged(s: &Summary) -> HashMap<usize, f64> {
1002 const K: usize = 6;
1004 let mut last: HashMap<u32, usize> = HashMap::new();
1005 for (i, c) in s.checks.iter().enumerate() {
1006 last.insert(c.tick, i);
1007 }
1008 let mut rows: Vec<(usize, Obs<K, 0>)> = Vec::new();
1009 for (ci, c) in s.checks.iter().enumerate() {
1010 for i in c.first as usize..(c.first as usize + c.n as usize).min(s.applies.len()) {
1011 let a = &s.applies[i];
1012 if model_route(a).is_none() || a.label == Label::Open {
1013 continue;
1014 }
1015 let y = f64::from(u8::from(matches!(a.label, Label::Changed | Label::Ignored)));
1016 let x = [
1017 1.0,
1018 f64::from(u8::from(a.route == Route::Ask)),
1019 f64::from(u8::from(a.route == Route::Free)),
1020 f64::from(u8::from(c.stop)),
1021 (1.0 + f64::from(a.pending)).ln(),
1022 f64::from(u8::from(last.get(&c.tick) == Some(&ci))),
1023 ];
1024 rows.push((i, Obs { x, y, w: 1.0, lv: [] }));
1025 }
1026 }
1027 if rows.is_empty() {
1028 return HashMap::new();
1029 }
1030 let spec = [Coef { mean: 0.0, sd: PRIOR_SD, ard: false, nonpos: false }; K];
1031 let obs: Vec<Obs<K, 0>> = rows.iter().map(|(_, o)| o.clone()).collect();
1032 let Some(est) = map(&obs, [], &spec, &[PRIOR_SD; K], &[], None) else { return HashMap::new() };
1033 rows.iter().map(|(i, o)| (*i, sigmoid(est.eta(o)))).collect()
1034}
1035
1036fn loglik(est: &Est<P, FAMILIES>, obs: &[Obs<P, FAMILIES>]) -> f64 {
1038 obs.iter()
1039 .map(|o| {
1040 let p = sigmoid(est.eta(o)).clamp(1e-12, 1.0 - 1e-12);
1042 o.w * (o.y * p.ln() + (1.0 - o.y) * (1.0 - p).ln())
1043 })
1044 .sum()
1045}
1046
1047fn half_life(d: &Data, base: &Est<P, FAMILIES>, spec: &[Coef; P]) -> Option<u32> {
1051 let mut ticks: Vec<u32> = d.tick.clone();
1052 ticks.sort_unstable();
1053 ticks.dedup();
1054 let eval: Vec<u32> = ticks.iter().skip(1).rev().take(PREQ_SESSIONS).rev().copied().collect();
1055 if eval.is_empty() {
1056 return None;
1057 }
1058 let mut best: (f64, Option<u32>) = (f64::NEG_INFINITY, None);
1059 for h in HALF_LIVES {
1060 let mut score = 0.0;
1061 let mut warm: Option<Est<P, FAMILIES>> = None;
1062 for &t in &eval {
1063 let train: Vec<Obs<P, FAMILIES>> = d
1064 .obs
1065 .iter()
1066 .zip(&d.tick)
1067 .zip(&d.base_w)
1068 .filter(|((_, tk), _)| **tk < t)
1069 .map(|((o, tk), w)| Obs { w: w * decay(t - tk, h), ..o.clone() })
1070 .collect();
1071 let test: Vec<Obs<P, FAMILIES>> = d.obs.iter().zip(&d.tick).filter(|(_, tk)| **tk == t).map(|(o, _)| o.clone()).collect();
1072 let Some(est) = map(&train, d.levels, spec, &base.sd, &base.spread, warm.as_ref()) else { continue };
1073 score += loglik(&est, &test);
1074 warm = Some(est);
1075 }
1076 if score > best.0 + 1e-9 * best.0.abs().max(1.0) || best.0 == f64::NEG_INFINITY {
1078 best = (score, h);
1079 }
1080 }
1081 best.1
1082}
1083
1084fn model_of(est: &Est<P, FAMILIES>, d: &Data, judge: String) -> Option<Model> {
1087 let mut chol: Vec<f64> = est.cov.iter().flatten().copied().collect();
1088 let finite = est.beta.iter().chain(&est.sd).chain(&est.spread).chain(est.u.iter().flatten()).chain(est.var_u.iter().flatten()).all(|x| x.is_finite());
1089 if !finite || !cholesky(&mut chol, P) {
1090 return None;
1091 }
1092 let levels = |f: usize| -> BTreeMap<String, (f64, f64)> {
1093 d.names[f].iter().enumerate().map(|(l, n)| (n.clone(), (est.u[f][l], est.var_u[f][l].max(0.0).sqrt()))).collect()
1094 };
1095 Some(Model {
1096 beta: est.beta.to_vec(),
1097 chol,
1098 sd: est.sd.to_vec(),
1099 spread: est.spread,
1100 sessions: levels(SESSION),
1101 lessons: levels(LESSON),
1102 judges: levels(JUDGE),
1103 judge,
1104 converged: est.converged,
1105 })
1106}
1107
1108#[must_use]
1111pub fn fit(s: &Summary, budget: Option<f64>) -> Fit {
1112 fit_with(s, budget, &coefs(), true)
1113}
1114
1115#[must_use]
1121pub fn fit_pruned(s: &Summary, budget: Option<f64>) -> Fit {
1122 fit_with(s, budget, &pruned(), false)
1123}
1124
1125fn fit_with(s: &Summary, budget: Option<f64>, spec: &[Coef; P], levels: bool) -> Fit {
1127 let data = |zero: bool| if levels { Data::of(s, zero) } else { Data::of(s, zero).without_levels() };
1128 let d = data(false);
1129 let now = s.sessions.len() as u32;
1130 let mut fit = Fit { key: s.key(), labels: d.obs.len() as u32, dropped: d.dropped, budget: budget.filter(|b| *b > 0.0), ..Fit::default() };
1131 let base = if d.obs.is_empty() { None } else { eb(&d.obs, d.levels, spec, None) };
1132 if let Some(base) = &base {
1133 fit.half_life = half_life(&d, base, spec);
1134 let est = eb(&d.weighted(now, fit.half_life), d.levels, spec, Some(base));
1135 let judge = s.checks.iter().rev().flat_map(|c| applies_of(s, c)).find_map(|a| a.judge);
1136 let judge = judge.filter(|_| levels).and_then(|j| s.judges.get(j as usize)).cloned().unwrap_or_default();
1137 if let Some(est) = est {
1138 let z = data(true);
1139 if z.unclear.iter().any(|u| *u)
1140 && let Some(alt) = eb(&z.weighted(now, fit.half_life), z.levels, spec, None)
1141 {
1142 fit.sensitive = [ASK, NOT].iter().any(|&k| est.beta[k].signum() != alt.beta[k].signum());
1143 }
1144 let eta = d.obs.iter().map(|o| est.eta(o)).sum::<f64>() / d.obs.len() as f64;
1145 fit.pool.eta = eta;
1146 fit.pool.slope = (est.beta[POOL], est.cov[POOL][POOL].max(0.0).sqrt());
1147 fit.pool.not = (est.beta[NOT], est.cov[NOT][NOT].max(0.0).sqrt());
1148 fit.model = model_of(&est, &d, judge).map(|m| if levels { m } else { Model { spread: [0.0; FAMILIES], ..m } });
1150 }
1151 }
1152 fit.rates = rates(s, fit.half_life);
1153 pool_state(s, &mut fit.pool);
1154 if let Some(b) = fit.budget {
1155 fit.lambda = price(s, &fit, b);
1156 }
1157 fit
1158}
1159
1160#[must_use]
1165pub fn prior() -> Fit {
1166 prior_with(&coefs(), true)
1167}
1168
1169#[must_use]
1173pub fn prior_pruned() -> Fit {
1174 prior_with(&pruned(), false)
1175}
1176
1177fn prior_with(spec: &[Coef; P], levels: bool) -> Fit {
1180 let s = Summary::default();
1181 let d = Data::of(&s, false);
1182 let sd: [f64; P] = std::array::from_fn(|k| spec[k].sd);
1183 let model = map(&d.obs, d.levels, spec, &sd, &[LEVEL_START; FAMILIES], None)
1184 .and_then(|est| model_of(&est, &d, String::new()))
1185 .map(|m| if levels { m } else { Model { spread: [0.0; FAMILIES], ..m } });
1186 Fit { model, rates: rates(&s, None), ..Fit::default() }
1187}
1188
1189const PROBIT: f64 = std::f64::consts::PI / 8.0;
1192
1193#[must_use]
1199pub fn predict(fit: &Fit, f: &Feat, route: Route, m: usize, session: &str, lesson: &str) -> Option<f64> {
1200 let model = fit.model.as_ref().filter(|m| m.current())?;
1201 let x = design(f, route, m);
1202 let lx = (0..P).map(|k| (k..P).map(|i| model.chol[i * P + k] * x[i]).sum::<f64>());
1204 let mut var: f64 = lx.map(|v| v * v).sum();
1205 let mut mean: f64 = model.beta.iter().zip(&x).map(|(b, x)| b * x).sum();
1206 for (map, name, spread) in [
1207 (&model.sessions, session, model.spread[SESSION]),
1208 (&model.lessons, lesson, model.spread[LESSON]),
1209 (&model.judges, model.judge.as_str(), model.spread[JUDGE]),
1210 ] {
1211 let (mu, sd) = map.get(name).copied().unwrap_or((0.0, spread));
1212 mean += mu;
1213 var += sd * sd;
1214 }
1215 Some(sigmoid(mean / (1.0 + PROBIT * var).sqrt()))
1216}
1217
1218fn rates(s: &Summary, h: Option<u32>) -> Rates {
1221 let mean = |(k, miss): (f64, f64)| (1.0 + k) / (2.0 + k + miss);
1222 let [free_stop, free_live, later, wait, ask] = s.tallies(h);
1223 let ask_live = if ask.0 + ask.1 > 0.0 { mean(ask) } else { 1.0 };
1225 Rates { free_stop: mean(free_stop), free_live: mean(free_live), later: mean(later), wait: mean(wait), ask_live }
1226}
1227
1228fn pool_state(s: &Summary, pool: &mut PoolFit) {
1231 let bins: Vec<(f64, f64)> = s.pool.bins.iter().take_while(|(n, _)| *n > 0).map(|&(n, k)| (1.0 + k as f64, 1.0 + n.saturating_sub(k) as f64)).collect();
1235 pool.yields = pava(&bins);
1236 pool.per_lesson = if s.pool.shown > 0 { s.pool.tokens as f64 / s.pool.shown as f64 } else { 0.0 };
1237 let obs = &s.pool.latency;
1238 pool.latency_n = obs.len() as u32;
1239 let mut read: Vec<u32> = obs.iter().map(|o| o.0).collect();
1240 read.sort_unstable_by(|a, b| b.cmp(a));
1241 pool.latency_reach = read.get(LATENCY_MIN - 1).copied().unwrap_or(0);
1242 if obs.len() >= 2 {
1243 let n = obs.len() as f64;
1244 let mx = obs.iter().map(|o| f64::from(o.0)).sum::<f64>() / n;
1245 let my = obs.iter().map(|o| f64::from(o.1)).sum::<f64>() / n;
1246 let sxx: f64 = obs.iter().map(|o| (f64::from(o.0) - mx).powi(2)).sum();
1247 let sxy: f64 = obs.iter().map(|o| (f64::from(o.0) - mx) * (f64::from(o.1) - my)).sum();
1248 let b = if sxx > 0.0 { (sxy / sxx).max(0.0) } else { 0.0 };
1249 let a = my - b * mx;
1250 let mut res: Vec<f64> = obs.iter().map(|o| f64::from(o.1) - a - b * f64::from(o.0)).collect();
1251 res.sort_by(f64::total_cmp);
1252 let q = res[((1.0 - DELTA) * (res.len() - 1) as f64).ceil() as usize];
1253 pool.latency = (a, b, q.max(0.0));
1254 }
1255}
1256
1257fn pava(bins: &[(f64, f64)]) -> Vec<(f64, f64)> {
1260 let mut blocks: Vec<(f64, f64, usize)> = Vec::new();
1261 for &(a, b) in bins {
1262 blocks.push((a, b, 1));
1263 while blocks.len() >= 2 {
1264 let (x, y) = (blocks[blocks.len() - 2], blocks[blocks.len() - 1]);
1265 if x.0 / (x.0 + x.1) >= y.0 / (y.0 + y.1) {
1266 break;
1267 }
1268 blocks.pop();
1269 if let Some(last) = blocks.last_mut() {
1270 *last = (x.0 + y.0, x.1 + y.1, x.2 + y.2);
1271 }
1272 }
1273 }
1274 blocks.iter().flat_map(|&(a, b, n)| std::iter::repeat_n((a, b), n)).collect()
1275}
1276
1277#[derive(Debug, Clone)]
1280struct Theta {
1281 b: [f64; P],
1282 base: Vec<f64>,
1283}
1284
1285struct Levels<'a> {
1288 name: &'a str,
1289 session: u64,
1290 judge: u64,
1291 lessons: Vec<u64>,
1292}
1293
1294impl<'a> Levels<'a> {
1295 fn of(session: &'a str, judge: &str, units: &[Unit]) -> Levels<'a> {
1296 Levels {
1297 name: session,
1298 session: hash(&[b"session", session.as_bytes()]),
1299 judge: hash(&[b"judge", judge.as_bytes()]),
1300 lessons: units.iter().map(|u| hash(&[b"lesson", u.lesson.as_bytes()])).collect(),
1301 }
1302 }
1303}
1304
1305fn draw_beta(m: &Model, seed: u64) -> [f64; P] {
1312 let mut rng = Rng::new(seed);
1313 let mut z: [f64; P] = std::array::from_fn(|_| normal(&mut rng));
1314 let spec = coefs();
1315 let l = |i: usize, k: usize| m.chol.get(i * P + k).copied().unwrap_or(0.0);
1316 let mut x = [0.0; P];
1317 for i in 0..P {
1318 let given = m.beta.get(i).copied().unwrap_or(0.0) + (0..i).map(|k| l(i, k) * z[k]).sum::<f64>();
1319 x[i] = given + l(i, i) * z[i];
1320 if spec[i].nonpos {
1321 x[i] = below_zero(given, l(i, i), z[i]);
1322 if l(i, i) > 0.0 {
1324 z[i] = (x[i] - given) / l(i, i);
1325 }
1326 }
1327 }
1328 x
1329}
1330
1331fn below_zero(mean: f64, sd: f64, z: f64) -> f64 {
1335 if sd <= 0.0 || !sd.is_finite() {
1336 return mean.min(0.0);
1337 }
1338 let bound = -mean / sd;
1339 let ln_lower = ln_phi(z) + ln_phi(bound);
1341 let upper = phi(-z) + phi(z) * phi(-bound);
1342 (mean + sd * probit(ln_lower, upper)).min(0.0)
1343}
1344
1345fn phi(x: f64) -> f64 {
1347 0.5 * erfc(-x / std::f64::consts::SQRT_2)
1348}
1349
1350fn ln_phi(x: f64) -> f64 {
1352 if x < 0.0 { 0.5f64.ln() - 0.5 * x * x + erfcx(-x / std::f64::consts::SQRT_2).ln() } else { (-0.5 * erfc(x / std::f64::consts::SQRT_2)).ln_1p() }
1353}
1354
1355fn erfc(x: f64) -> f64 {
1357 let a = x.abs();
1358 let tail = (-a * a).exp() * erfcx(a);
1359 if x >= 0.0 { tail } else { 2.0 - tail }
1360}
1361
1362fn erfcx(x: f64) -> f64 {
1365 const C: [f64; 10] =
1366 [-1.265_512_23, 1.000_023_68, 0.374_091_96, 0.096_784_18, -0.186_288_06, 0.278_868_07, -1.135_203_98, 1.488_515_87, -0.822_152_23, 0.170_872_77];
1367 let t = 1.0 / (1.0 + 0.5 * x);
1368 t * C.iter().rev().fold(0.0, |acc, c| acc * t + c).exp()
1369}
1370
1371fn probit(ln_lower: f64, upper: f64) -> f64 {
1375 const A: [f64; 6] =
1376 [-3.969_683_028_665_376e1, 2.209_460_984_245_205e2, -2.759_285_104_469_687e2, 1.383_577_518_672_69e2, -3.066_479_806_614_716e1, 2.506_628_277_459_239];
1377 const B: [f64; 5] = [-5.447_609_879_822_406e1, 1.615_858_368_580_409e2, -1.556_989_798_598_866e2, 6.680_131_188_771_972e1, -1.328_068_155_288_572e1];
1378 const C: [f64; 6] =
1379 [-7.784_894_002_430_293e-3, -3.223_964_580_411_365e-1, -2.400_758_277_161_838, -2.549_732_539_343_734, 4.374_664_141_464_968, 2.938_163_982_698_783];
1380 const D: [f64; 4] = [7.784_695_709_041_462e-3, 3.224_671_290_700_398e-1, 2.445_134_137_142_996, 3.754_408_661_907_416];
1381 const LOW: f64 = 0.02425;
1382 let poly = |c: &[f64], x: f64| c.iter().fold(0.0, |acc, k| acc * x + k);
1383 let tail = |ln_p: f64| {
1385 let q = (-2.0 * ln_p).sqrt();
1386 poly(&C, q) / (poly(&D, q) * q + 1.0)
1387 };
1388 let lower = ln_lower.exp();
1389 if lower < LOW {
1390 tail(ln_lower)
1391 } else if upper < LOW {
1392 -tail(upper.ln())
1393 } else {
1394 let q = lower - 0.5;
1395 let r = q * q;
1396 poly(&A, r) * q / (poly(&B, r) * r + 1.0)
1397 }
1398}
1399
1400fn theta_with(m: &Model, b: [f64; P], feats: &[Feat], units: &[Unit], lv: &Levels, seed: u64) -> Theta {
1405 theta_by(m, b, feats, units, lv, |(mean, sd), h| mean + sd * normal(&mut Rng::new(mix(seed, h))))
1406}
1407
1408fn theta_mean(m: &Model, feats: &[Feat], units: &[Unit], lv: &Levels) -> Theta {
1414 let spec = coefs();
1415 let b = std::array::from_fn(|i| {
1416 let x = m.beta.get(i).copied().unwrap_or(0.0);
1417 if spec[i].nonpos { x.min(0.0) } else { x }
1418 });
1419 theta_by(m, b, feats, units, lv, |(mean, _), _| mean)
1420}
1421
1422fn theta_by(m: &Model, b: [f64; P], feats: &[Feat], units: &[Unit], lv: &Levels, level: impl Fn((f64, f64), u64) -> f64) -> Theta {
1425 let effect = |map: &BTreeMap<String, (f64, f64)>, name: &str, spread: f64, h: u64| level(map.get(name).copied().unwrap_or((0.0, spread)), h);
1426 let shared = effect(&m.sessions, lv.name, m.spread[SESSION], lv.session) + effect(&m.judges, &m.judge, m.spread[JUDGE], lv.judge);
1427 let base = feats
1428 .iter()
1429 .zip(units)
1430 .zip(&lv.lessons)
1431 .map(|((f, u), h)| dot(&b, &design(f, Route::Free, 1)) + shared + effect(&m.lessons, u.lesson, m.spread[LESSON], *h))
1432 .collect();
1433 Theta { b, base }
1434}
1435
1436fn theta(m: &Model, feats: &[Feat], units: &[Unit], lv: &Levels, seed: u64) -> Theta {
1438 theta_with(m, draw_beta(m, seed), feats, units, lv, seed)
1439}
1440
1441struct Valuer<'a> {
1443 t: &'a Theta,
1444 feats: &'a [Feat],
1445 rates: &'a Rates,
1446}
1447
1448impl Valuer<'_> {
1449 fn n(&self) -> usize {
1450 self.feats.len()
1451 }
1452
1453 fn stop(&self) -> bool {
1454 self.feats.first().is_some_and(|f| f.stop)
1455 }
1456
1457 fn p_not(&self, j: usize) -> f64 {
1459 let b = &self.t.b;
1460 sigmoid(self.t.base[j] + b[NOT] + b[NOT_G] * (f64::from(self.feats[j].grade) - 2.0))
1461 }
1462
1463 fn ask(&self, j: usize, m: usize) -> f64 {
1464 let (b, f) = (&self.t.b, &self.feats[j]);
1465 let eta = self.t.base[j]
1466 + b[ASK]
1467 + b[ASK_M] * (m.max(1) as f64).ln()
1468 + b[ASK_STOP] * f64::from(u8::from(f.stop))
1469 + b[ASK_G3] * f64::from(u8::from(f.grade == 3));
1470 let rate = match (self.stop(), f.three) {
1475 (true, _) => 1.0,
1476 (false, true) => self.rates.ask_live,
1477 (false, false) => self.rates.free_live,
1478 };
1479 rate * (sigmoid(eta) - self.p_not(j))
1480 }
1481
1482 fn free(&self, j: usize, m: usize) -> f64 {
1483 let rate = if self.stop() { self.rates.free_stop } else { self.rates.free_live };
1484 rate * (sigmoid(self.t.base[j] + self.t.b[FREE_M] * (m.max(1) as f64).ln()) - self.p_not(j))
1485 }
1486
1487 fn value(&self, routes: &[Route]) -> f64 {
1493 let m_q = routes.iter().filter(|r| **r == Route::Ask).count();
1494 let m_f = routes.iter().filter(|r| **r == Route::Free).count();
1495 routes
1496 .iter()
1497 .enumerate()
1498 .map(|(j, r)| match r {
1499 Route::Ask => self.ask(j, m_q),
1500 Route::Free => self.free(j, m_f),
1501 Route::Placebo | Route::Not => 0.0,
1502 })
1503 .sum()
1504 }
1505}
1506
1507fn forced(f: &Feat) -> bool {
1510 f.stop && f.grade == 3
1511}
1512
1513struct Table {
1516 n: usize,
1517 ask: Vec<f64>,
1521 free: Vec<f64>,
1522 forced: Vec<bool>,
1523 nofree: Vec<bool>,
1524}
1525
1526impl Table {
1527 fn of(v: &Valuer, nofree: &[bool]) -> Table {
1530 let n = v.n();
1531 let mut ask = vec![0.0; (n + 1) * n];
1532 let mut free = vec![0.0; (n + 1) * n];
1533 for m in 1..=n {
1534 for j in 0..n {
1535 ask[m * n + j] = v.ask(j, m);
1536 free[m * n + j] = v.free(j, m);
1537 }
1538 }
1539 Table { n, ask, free, forced: v.feats.iter().map(forced).collect(), nofree: (0..n).map(|j| nofree.get(j).copied().unwrap_or(false)).collect() }
1540 }
1541
1542 #[cfg(test)]
1545 fn value(&self, routes: &[Route]) -> f64 {
1546 let m_q = routes.iter().filter(|r| **r == Route::Ask).count();
1547 let m_f = routes.iter().filter(|r| **r == Route::Free).count();
1548 let n = self.n;
1549 routes
1550 .iter()
1551 .enumerate()
1552 .map(|(j, r)| match r {
1553 Route::Ask => self.ask[m_q * n + j],
1554 Route::Free => self.free[m_f * n + j],
1555 Route::Placebo | Route::Not => 0.0,
1556 })
1557 .sum()
1558 }
1559}
1560
1561#[derive(Debug, Clone, Copy, PartialEq)]
1563struct Total(f64);
1564
1565impl Eq for Total {}
1566
1567impl PartialOrd for Total {
1568 fn partial_cmp(&self, other: &Self) -> Option<std::cmp::Ordering> {
1569 Some(self.cmp(other))
1570 }
1571}
1572
1573impl Ord for Total {
1574 fn cmp(&self, other: &Self) -> std::cmp::Ordering {
1575 self.0.total_cmp(&other.0)
1576 }
1577}
1578
1579fn top_sums(values: impl Iterator<Item = f64>, k: usize, out: &mut Vec<f64>, heap: &mut std::collections::BinaryHeap<std::cmp::Reverse<Total>>) {
1582 out.clear();
1583 heap.clear();
1584 out.push(if k == 0 { 0.0 } else { f64::NEG_INFINITY });
1585 let mut sum = 0.0;
1586 for x in values {
1587 if k > 0 {
1588 if heap.len() < k {
1589 heap.push(std::cmp::Reverse(Total(x)));
1590 sum += x;
1591 } else if let Some(min) = heap.peek().map(|r| r.0.0)
1592 && x > min
1593 {
1594 heap.pop();
1595 heap.push(std::cmp::Reverse(Total(x)));
1596 sum += x - min;
1597 }
1598 }
1599 out.push(if heap.len() < k { f64::NEG_INFINITY } else { sum });
1600 }
1601}
1602
1603fn solve(t: &Table, rail: usize, cap: usize) -> Vec<Route> {
1611 let n = t.n;
1612 let nf_forced = t.forced.iter().filter(|x| **x).count();
1613 let open: Vec<usize> = (0..n).filter(|&j| !t.forced[j]).collect();
1614 let held = open.iter().filter(|&&j| t.nofree[j]).count();
1617 let q_max = rail.max(nf_forced).min(n);
1618 let f_max = cap.min(open.len() - held);
1619 let ask = |m: usize, j: usize| if m == 0 { 0.0 } else { t.ask[m * n + j] };
1621 let free = |m: usize, j: usize| if m == 0 { 0.0 } else { t.free[m * n + j] };
1622 let mut order = open.clone();
1623 let mut keys = vec![0.0; n];
1624 let (mut pre, mut suf) = (Vec::with_capacity(n + 1), Vec::with_capacity(n + 1));
1625 let mut heap = std::collections::BinaryHeap::with_capacity(n + 1);
1626 let sort = |order: &mut Vec<usize>, keys: &mut Vec<f64>, mq: usize, mf: usize| {
1627 for &j in &open {
1628 keys[j] = if t.nofree[j] { f64::INFINITY } else { ask(mq, j) - free(mf, j) };
1629 }
1630 order.clear();
1631 order.extend_from_slice(&open);
1632 order.sort_by(|&a, &b| keys[b].total_cmp(&keys[a]).then(a.cmp(&b)));
1633 };
1634 let mut best: Option<(f64, usize, usize, usize)> = None;
1635 for mq in nf_forced..=q_max {
1636 let k = mq - nf_forced;
1637 if k > open.len() {
1638 break;
1639 }
1640 let fixed: f64 = (0..n).filter(|&j| t.forced[j]).map(|j| ask(mq, j)).sum();
1641 for mf in 0..=f_max.min(open.len() - k) {
1642 sort(&mut order, &mut keys, mq, mf);
1643 top_sums(order.iter().map(|&j| ask(mq, j)), k, &mut pre, &mut heap);
1644 top_sums(order.iter().rev().map(|&j| free(mf, j)), mf, &mut suf, &mut heap);
1645 for split in k.max(held)..=open.len() - mf {
1646 let total = fixed + pre[split] + suf[open.len() - split];
1647 if best.is_none_or(|b| total > b.0) {
1648 best = Some((total, mq, mf, split));
1649 }
1650 }
1651 }
1652 }
1653 let mut routes = vec![Route::Not; n];
1654 for j in (0..n).filter(|&j| t.forced[j]) {
1655 routes[j] = Route::Ask;
1656 }
1657 let Some((_, mq, mf, split)) = best else { return routes };
1658 sort(&mut order, &mut keys, mq, mf);
1659 let pick = |within: &[usize], gain: &dyn Fn(usize) -> f64, k: usize| -> Vec<usize> {
1660 let mut c = within.to_vec();
1661 c.sort_by(|&a, &b| gain(b).total_cmp(&gain(a)).then(a.cmp(&b)));
1662 c.truncate(k);
1663 c
1664 };
1665 for j in pick(&order[..split], &|j| ask(mq, j), mq - nf_forced) {
1666 routes[j] = Route::Ask;
1667 }
1668 for j in pick(&order[split..], &|j| free(mf, j), mf) {
1669 routes[j] = Route::Free;
1670 }
1671 routes
1672}
1673
1674#[must_use]
1677pub fn fixed_routes(units: &[Unit], had: &[String]) -> Vec<Route> {
1678 let ls: Vec<Challenge> = units.iter().map(|u| Challenge { id: u.lesson.to_string(), ..Default::default() }).collect();
1679 let applied: Vec<(&Challenge, u8)> = ls.iter().zip(units).map(|(c, u)| (c, u.grade)).collect();
1680 let (ask, free) = conscience::graded_without(&applied, None, had);
1681 let holds = |v: &[&Challenge], c: &Challenge| v.iter().any(|x| std::ptr::eq(*x, c));
1682 ls.iter()
1683 .map(|c| match (holds(&ask, c), holds(&free, c)) {
1684 (true, _) => Route::Ask,
1685 (_, true) => Route::Free,
1686 _ => Route::Not,
1687 })
1688 .collect()
1689}
1690
1691fn nofree_of(units: &[Unit], had: &[String]) -> Vec<bool> {
1693 units.iter().map(|u| had.iter().any(|h| h == u.lesson)).collect()
1694}
1695
1696#[derive(Debug, Clone, Copy)]
1698pub struct Check<'a> {
1699 pub session: &'a str,
1700 pub stop: bool,
1701 pub tokens: u32,
1703}
1704
1705pub const FIXED: &str = "fixed";
1709pub const SIZE: &str = "size";
1711pub const GATE: &str = "gate";
1714pub const EVIDENCE: &str = "evidence";
1718pub const EXPERIMENT: &str = "experiment";
1721pub const EXPERIMENT_CONTROL: &str = "experiment_control";
1724
1725#[derive(Debug, Clone, PartialEq)]
1736pub struct Decision {
1737 pub routes: Vec<Route>,
1738 pub p: Vec<f64>,
1739 pub pf: Vec<f64>,
1740 pub source: &'static str,
1741 pub gate: f64,
1742 pub risk: f64,
1743 pub xk: Option<u32>,
1744 pub xs: Option<u32>,
1745}
1746
1747fn lower_quantile(v: &mut [f64]) -> f64 {
1749 if v.is_empty() {
1750 return 0.0;
1751 }
1752 v.sort_by(f64::total_cmp);
1753 v[((DELTA * v.len() as f64).ceil() as usize).clamp(1, v.len()) - 1]
1754}
1755
1756#[cfg(test)]
1759#[must_use]
1760pub(crate) fn remaining(gate: &Gate) -> f64 {
1761 gate.allowance + lower_quantile(&mut gate.gain.clone())
1762}
1763
1764fn allowance_at(gate: &Gate, alpha: f64) -> f64 {
1767 if alpha > 0.0 { gate.allowance * (alpha / ALPHA) } else { 0.0 }
1770}
1771
1772#[must_use]
1775pub fn remaining_at(gate: &Gate, alpha: f64) -> f64 {
1776 let (allowance, ledger) = gate_parts(gate, alpha);
1777 (allowance + ledger).max(0.0)
1778}
1779
1780#[must_use]
1785pub fn gate_parts(gate: &Gate, alpha: f64) -> (f64, f64) {
1786 (allowance_at(gate, alpha), lower_quantile(&mut gate.gain.clone()))
1787}
1788
1789fn gate_seed(fit: &Fit, k: usize) -> u64 {
1792 mix(hash(&[b"gate", fit.key.as_bytes()]), k as u64)
1793}
1794
1795struct Draws<'a> {
1798 thetas: Vec<Theta>,
1799 feats: Vec<Feat>,
1800 rates: &'a Rates,
1801 base: Vec<f64>,
1802}
1803
1804impl<'a> Draws<'a> {
1805 fn of(model: &Model, fit: &'a Fit, units: &[Unit], check: &Check, lv: &Levels, base: &[Route]) -> Draws<'a> {
1806 let feats = Feat::of(units, check.stop, check.tokens);
1807 let thetas: Vec<Theta> = (0..DRAWS).map(|k| theta(model, &feats, units, lv, gate_seed(fit, k))).collect();
1808 let mut d = Draws { thetas, feats, rates: &fit.rates, base: Vec::new() };
1809 d.base = (0..DRAWS).map(|k| d.value(k, base)).collect();
1810 d
1811 }
1812
1813 fn value(&self, k: usize, routes: &[Route]) -> f64 {
1815 Valuer { t: &self.thetas[k], feats: &self.feats, rates: self.rates }.value(routes)
1816 }
1817
1818 fn gains(&self, routes: &[Route]) -> Vec<f64> {
1820 self.base.iter().enumerate().map(|(k, b)| self.value(k, routes) - b).collect()
1821 }
1822}
1823
1824fn evident(model: &Model, d: &Draws, units: &[Unit], lv: &Levels, base: &[Route], nofree: &[bool]) -> Option<Vec<Route>> {
1828 let t = theta_mean(model, &d.feats, units, lv);
1829 let best = solve(&Table::of(&Valuer { t: &t, feats: &d.feats, rates: d.rates }, nofree), crate::attention::RAIL, conscience::CAP);
1830 (best != base && lower_quantile(&mut d.gains(&best)) > 0.0).then_some(best)
1831}
1832
1833pub(crate) fn tested(units: &[Unit], base: &[Route], had: &[String]) -> Vec<usize> {
1838 let mut order: Vec<usize> =
1839 (0..units.len()).filter(|&j| units[j].grade != 3 && matches!(base[j], Route::Ask | Route::Free) && !had.iter().any(|h| h == units[j].lesson)).collect();
1840 order.sort_by_key(|&j| std::cmp::Reverse(units[j].grade));
1841 order
1842}
1843
1844fn test_routes(base: &[Route], order: &[usize], s: usize) -> Vec<Route> {
1847 let mut routes = base.to_vec();
1848 for (i, &j) in order.iter().enumerate() {
1849 routes[j] = if i < s { Route::Ask } else { Route::Free };
1850 }
1851 routes
1852}
1853
1854fn draw_count(seed: u64, k: usize) -> usize {
1858 (hash(&[b"question-count", &seed.to_le_bytes()]) % (k as u64 + 1)) as usize
1859}
1860
1861#[must_use]
1877pub fn choose(units: &[Unit], check: &Check, fit: &Fit, gate: &Gate, alpha: f64, seed: u64, had: &[String]) -> Decision {
1878 let base = fixed_routes(units, had);
1879 let fixed = |source: &'static str| Decision {
1880 p: share(&[&base], Route::Ask),
1881 pf: share(&[&base], Route::Free),
1882 routes: base.clone(),
1883 source,
1884 gate: 0.0,
1885 risk: 0.0,
1886 xk: None,
1887 xs: None,
1888 };
1889 let Some(model) = fit.model.as_ref().filter(|m| m.current()) else { return fixed(FIXED) };
1890 if gate.allowance <= 0.0 || gate.gain.len() != DRAWS || units.is_empty() {
1891 return fixed(FIXED);
1892 }
1893 if units.len() > SOLVE_MAX {
1894 return fixed(SIZE);
1895 }
1896 let lv = Levels::of(check.session, &model.judge, units);
1897 let d = Draws::of(model, fit, units, check, &lv, &base);
1898 let allowance = allowance_at(gate, alpha);
1899 let held = lower_quantile(&mut gate.gain.clone());
1900 let remaining = remaining_at(gate, alpha);
1901 let admit = |routes: &[Route]| -> Option<f64> {
1905 if routes == base {
1906 return Some(0.0);
1907 }
1908 let mut after: Vec<f64> = d.gains(routes).iter().zip(&gate.gain).map(|(x, g)| g + x).collect();
1909 let q = lower_quantile(&mut after);
1910 (q >= -allowance).then_some((held - q).max(0.0))
1911 };
1912 if let Some(best) = evident(model, &d, units, &lv, &base, &nofree_of(units, had)) {
1913 return match admit(&best) {
1914 Some(risk) => Decision {
1915 p: share(&[&best], Route::Ask),
1916 pf: share(&[&best], Route::Free),
1917 routes: best,
1918 source: EVIDENCE,
1919 gate: remaining,
1920 risk,
1921 xk: None,
1922 xs: None,
1923 },
1924 None => Decision { gate: remaining, ..fixed(GATE) },
1925 };
1926 }
1927 let order = tested(units, &base, had);
1928 if check.stop || alpha <= 0.0 || order.is_empty() {
1929 return Decision { gate: remaining, ..fixed(EVIDENCE) };
1930 }
1931 let k = order.len();
1934 let arms: Vec<(Vec<Route>, Option<f64>)> = (0..=k)
1935 .map(|s| {
1936 let routes = test_routes(&base, &order, s);
1937 let risk = admit(&routes);
1938 (if risk.is_some() { routes } else { base.clone() }, risk)
1939 })
1940 .collect();
1941 let played: Vec<&Vec<Route>> = arms.iter().map(|(r, _)| r).collect();
1942 let s = draw_count(seed, k);
1943 let (routes, risk) = &arms[s];
1944 let source = match risk {
1945 None => GATE,
1946 Some(_) if *routes == base => EXPERIMENT_CONTROL,
1947 Some(_) => EXPERIMENT,
1948 };
1949 let count = |x: usize| Some(u32::try_from(x).unwrap_or(u32::MAX));
1950 Decision {
1951 p: share(&played, Route::Ask),
1952 pf: share(&played, Route::Free),
1953 routes: routes.clone(),
1954 source,
1955 gate: remaining,
1956 risk: risk.unwrap_or(0.0),
1957 xk: count(k),
1958 xs: count(s),
1959 }
1960}
1961
1962#[must_use]
1971pub fn decide_on_evidence(units: &[Unit], check: &Check, fit: &Fit, had: &[String]) -> Vec<Route> {
1972 let base = fixed_routes(units, had);
1973 let Some(model) = fit.model.as_ref().filter(|m| m.current()) else { return base };
1974 if units.is_empty() || units.len() > SOLVE_MAX {
1975 return base;
1976 }
1977 let lv = Levels::of(check.session, &model.judge, units);
1978 let d = Draws::of(model, fit, units, check, &lv, &base);
1979 evident(model, &d, units, &lv, &base, &nofree_of(units, had)).unwrap_or(base)
1980}
1981
1982fn share(draws: &[&Vec<Route>], route: Route) -> Vec<f64> {
1985 let n = draws.first().map_or(0, |d| d.len());
1986 (0..n).map(|j| crate::outcome::round_p(draws.iter().filter(|d| d[j] == route).count() as f64 / draws.len() as f64)).collect()
1987}
1988
1989pub fn ledger_add_without(fit: &Fit, gate: &mut Gate, check: &Check, units: &[Unit], routes: &[Route], had: &[String]) {
1995 let Some(model) = fit.model.as_ref().filter(|m| m.current()) else { return };
1996 let betas: Vec<[f64; P]> = (0..DRAWS).map(|k| draw_beta(model, gate_seed(fit, k))).collect();
1997 add_with(fit, &betas, gate, check, units, routes, had);
1998}
1999
2000fn add_with(fit: &Fit, betas: &[[f64; P]], gate: &mut Gate, check: &Check, units: &[Unit], routes: &[Route], had: &[String]) {
2003 let Some(model) = fit.model.as_ref().filter(|m| m.current()) else { return };
2004 if units.is_empty() || units.len() != routes.len() {
2005 return;
2006 }
2007 if gate.gain.len() != DRAWS {
2008 gate.gain = vec![0.0; DRAWS];
2009 }
2010 let base = fixed_routes(units, had);
2011 let routes: Vec<Route> = routes.iter().map(|r| if *r == Route::Placebo { Route::Not } else { *r }).collect();
2013 let feats = Feat::of(units, check.stop, check.tokens);
2014 let lv = Levels::of(check.session, &model.judge, units);
2015 let mut mean = 0.0;
2016 for (k, (g, b)) in gate.gain.iter_mut().zip(betas).enumerate() {
2017 let t = theta_with(model, *b, &feats, units, &lv, gate_seed(fit, k));
2018 let v = Valuer { t: &t, feats: &feats, rates: &fit.rates };
2019 let fixed = v.value(&base);
2020 *g += v.value(&routes) - fixed;
2021 mean += fixed / DRAWS as f64;
2022 }
2023 gate.allowance += ALPHA * mean.max(0.0);
2024 gate.checks += 1;
2025 if routes != base {
2026 gate.deviated += 1;
2027 }
2028}
2029
2030#[must_use]
2033pub fn ledger(s: &Summary) -> Gate {
2034 let mut gate = Gate::default();
2035 let Some(model) = s.fit.model.as_ref().filter(|m| m.current()) else { return gate };
2036 let betas: Vec<[f64; P]> = (0..DRAWS).map(|k| draw_beta(model, gate_seed(&s.fit, k))).collect();
2037 for c in s.checks.iter().filter(|c| c.learned) {
2038 let (units, _, routes) = check_view(s, c);
2039 let session = s.sessions.get(c.tick as usize).map_or("", |x| x.id.as_str());
2040 add_with(&s.fit, &betas, &mut gate, &Check { session, stop: c.stop, tokens: c.tokens }, &units, &routes, &logged_had(s, c));
2041 }
2042 gate
2043}
2044
2045#[derive(Debug, Clone, Default, PartialEq, Serialize, Deserialize)]
2047pub struct Monitor {
2048 #[serde(default, skip_serializing_if = "std::ops::Not::not")]
2050 pub reverted: bool,
2051 #[serde(default, skip_serializing_if = "is_default")]
2053 pub sessions: u32,
2054 #[serde(default, skip_serializing_if = "is_default")]
2055 pub units: u32,
2056 #[serde(default, skip_serializing_if = "is_default")]
2060 pub mean: f64,
2061 #[serde(default, skip_serializing_if = "Option::is_none")]
2062 pub ucb: Option<f64>,
2063 #[serde(default, skip_serializing_if = "is_default")]
2070 pub labelled: u32,
2071 #[serde(default, skip_serializing_if = "is_default")]
2072 pub ess: f64,
2073 #[serde(default, skip_serializing_if = "is_default")]
2074 pub max_weight: f64,
2075 #[serde(default, skip_serializing_if = "is_default")]
2076 pub truncated: u32,
2077}
2078
2079#[must_use]
2088pub fn monitor(s: &Summary) -> Monitor {
2089 let spec = coefs();
2090 let d = Data::of(s, false);
2091 let Some(model) = &s.fit.model else { return Monitor::default() };
2092 let sd: [f64; P] = std::array::from_fn(|k| model.sd.get(k).copied().unwrap_or(spec[k].sd));
2093 let spread = model.spread;
2094 let folds: Vec<Option<Est<P, FAMILIES>>> = (0..2u32)
2095 .map(|half| {
2096 let train: Vec<Obs<P, FAMILIES>> =
2097 d.weighted(s.sessions.len() as u32, s.fit.half_life).into_iter().zip(&d.tick).filter(|(_, t)| *t % 2 != half).map(|(o, _)| o).collect();
2098 map(&train, d.levels, &spec, &sd, &spread, None)
2099 })
2100 .collect();
2101 struct Term {
2103 tick: u32,
2104 diff: [f64; ARMS],
2105 p_hat: [f64; ARMS],
2106 taken: usize,
2107 y: Option<f64>,
2108 w: f64,
2109 vb: f64,
2110 }
2111 let mut terms: Vec<Term> = Vec::new();
2112 for c in s.checks.iter().filter(|c| c.learned) {
2113 let Some(est) = &folds[(c.tick % 2) as usize] else { continue };
2114 let (units, feats, routes) = check_view(s, c);
2115 let base = fixed_routes(&units, &logged_had(s, c));
2117 let (m_q, m_f) = (usize::from(c.m_q).max(1), usize::from(c.m_f).max(1));
2118 for (j, f) in feats.iter().enumerate() {
2119 let i = c.first as usize + j;
2120 let a = &s.applies[i];
2121 let lv = d.lv_of.get(&i).copied().unwrap_or([None; FAMILIES]);
2122 let eta = |route: Route, m: usize| -> f64 {
2124 dot(&est.beta, &design(f, route, m))
2125 + lv[1].map_or(0.0, |l| est.u[LESSON].get(l).copied().unwrap_or(0.0))
2126 + lv[2].map_or(0.0, |l| est.u[JUDGE].get(l).copied().unwrap_or(0.0))
2127 };
2128 let p_hat = [sigmoid(eta(Route::Ask, m_q)), sigmoid(eta(Route::Free, m_f)), sigmoid(eta(Route::Not, 1))];
2129 let mu = [a.p, a.pf, (1.0 - a.p - a.pf).max(0.0)];
2130 let (taken, fixed) = (routes[j].arm(), base[j].arm());
2131 let diff: [f64; ARMS] = std::array::from_fn(|x| mu[x] - f64::from(u8::from(x == fixed)));
2132 let y = match a.label {
2133 Label::Changed => Some(1.0),
2134 Label::Ignored => Some(0.0),
2135 _ => None,
2136 };
2137 let psi = if a.route == Route::Placebo { a.psi.max(PSI_FLOOR) } else { 1.0 };
2140 let judged_p = d.judged_p.get(&i).copied().unwrap_or(1.0).max(1.0 / MAX_IPW);
2141 let w = if y.is_some() { 1.0 / (mu[taken].max(1.0 / DRAWS as f64) * judged_p * psi) } else { 0.0 };
2142 let vb = p_hat[fixed] - p_hat[2];
2143 terms.push(Term { tick: c.tick, diff, p_hat, taken, y, w, vb });
2144 }
2145 }
2146 let labelled: Vec<f64> = terms.iter().filter(|t| t.y.is_some()).map(|t| t.w).collect();
2147 if terms.is_empty() || labelled.is_empty() {
2148 return Monitor::default();
2149 }
2150 let cap = (labelled.len() as f64).sqrt() * labelled.iter().sum::<f64>() / labelled.len() as f64;
2151 let mut norm = [0.0f64; ARMS];
2152 for t in &terms {
2153 norm[t.taken] += t.w.min(cap);
2154 }
2155 let count = terms.len() as f64;
2156 let norm = norm.map(|x| if x > 0.0 { x / count } else { 1.0 });
2157 let mut routed = [0usize; ARMS];
2161 for t in terms.iter().filter(|t| t.y.is_some() && t.diff[t.taken] != 0.0) {
2162 routed[t.taken] += 1;
2163 }
2164 let residual = routed.map(|n| n >= ROUTE_LABELS_MIN);
2165 let mut per: BTreeMap<u32, (f64, usize)> = BTreeMap::new();
2166 for t in &terms {
2167 let w = t.w.min(cap) / norm[t.taken];
2168 let r = |a: usize| if a == t.taken && residual[a] { t.y.map_or(0.0, |y| w * (y - t.p_hat[a])) } else { 0.0 };
2169 let dr: f64 = (0..ARMS).map(|a| t.diff[a] * (t.p_hat[a] + r(a))).sum();
2170 let e = per.entry(t.tick).or_insert((0.0, 0));
2171 e.0 += dr + ALPHA * t.vb;
2172 e.1 += 1;
2173 }
2174 let z: Vec<f64> = per.values().map(|(x, n)| x / *n as f64).collect();
2175 let mean = z.iter().sum::<f64>() / z.len() as f64;
2176 let ucb = Some(upper_bound(&z)).filter(|u| u.is_finite());
2177 let moves = |t: &&Term| t.y.is_some() && residual[t.taken] && t.diff[t.taken] != 0.0;
2181 let used: Vec<f64> = terms.iter().filter(moves).map(|t| t.w.min(cap) / norm[t.taken]).collect();
2182 let (sum, sq) = used.iter().fold((0.0, 0.0), |(s, q), w| (s + w, q + w * w));
2183 let ess = if sq > 0.0 { sum * sum / sq } else { 0.0 };
2184 let max_weight = if sum > 0.0 { used.iter().copied().fold(0.0, f64::max) * used.len() as f64 / sum } else { 0.0 };
2185 Monitor {
2186 reverted: ucb.is_some_and(|u| u < 0.0),
2187 sessions: z.len() as u32,
2188 units: terms.len() as u32,
2189 mean,
2190 ucb,
2191 labelled: used.len() as u32,
2192 ess,
2193 max_weight,
2194 truncated: terms.iter().filter(moves).filter(|t| t.w > cap).count() as u32,
2195 }
2196}
2197
2198fn upper_bound(x: &[f64]) -> f64 {
2207 let log = (1.0 / DELTA).ln();
2208 let shape = 2.0 * log + (1.0 + 2.0 * log).ln();
2209 let (mut sum, mut sum_sq) = (0.0, 0.0);
2210 let mut best = f64::INFINITY;
2211 for (i, &xi) in x.iter().enumerate() {
2212 let t = (i + 1) as f64;
2213 sum += xi;
2214 sum_sq += xi * xi;
2215 if i + 1 < CS_MIN.max(2) {
2216 continue;
2217 }
2218 let mean = sum / t;
2219 let var = ((sum_sq - t * mean * mean) / (t - 1.0)).max(f64::MIN_POSITIVE);
2220 let (v, rho) = (t * var, var * CS_SESSIONS / shape);
2221 best = best.min(mean + ((v + rho) * ((1.0 + v / rho).ln() + 2.0 * log)).sqrt() / t);
2222 }
2223 best
2224}
2225
2226#[must_use]
2228pub fn session_seed(session: &str) -> u64 {
2229 hash(&[b"pool", session.as_bytes()])
2230}
2231
2232#[must_use]
2234pub fn pool_grid() -> Vec<usize> {
2235 POOL_STEPS.iter().map(|&i| (crate::whole::TOKENS as f64 * 2f64.powf(f64::from(i) / 2.0)).round() as usize).collect()
2236}
2237
2238fn pool_value(pool: &PoolFit, tokens: usize, seed: u64) -> f64 {
2242 let mut rng = Rng::new(seed);
2243 let shown = tokens as f64 / pool.per_lesson.max(1.0);
2244 let mut applies = 0.0;
2245 for (b, &(a, bb)) in pool.yields.iter().enumerate() {
2246 let (lo, hi) = ((1u64 << b) - 1, (1u64 << (b + 1)) - 1);
2247 let inside = (shown.min(hi as f64) - lo as f64).max(0.0);
2248 applies += inside * if inside > 0.0 { beta_draw(a, bb, &mut rng) } else { 0.0 };
2249 }
2250 let slope = (pool.slope.0 + pool.slope.1 * normal(&mut rng)).min(0.0);
2251 let not = pool.not.0 + pool.not.1 * normal(&mut rng);
2252 let eta = pool.eta + slope * (tokens as f64 / crate::whole::TOKENS as f64).ln();
2253 applies * (sigmoid(eta) - sigmoid(eta + not)).max(0.0)
2254}
2255
2256fn pool_feasible(pool: &PoolFit, tokens: usize, limit_ms: Option<u64>, lambda: f64) -> bool {
2261 let baseline = crate::whole::TOKENS;
2262 let bins_to = |t: usize| -> usize {
2263 let shown = (t as f64 / pool.per_lesson.max(1.0)).ceil() as u64;
2264 (u64::BITS - shown.leading_zeros()) as usize
2265 };
2266 let known = pool.per_lesson > 0.0 && pool.yields.len() >= bins_to(tokens.max(baseline));
2267 if tokens > baseline {
2268 let (a, b, q) = pool.latency;
2269 let fits = limit_ms.is_some_and(|l| a + b * tokens as f64 + q <= l as f64);
2270 let observed = f64::from(pool.latency_reach) >= LATENCY_NEAR * tokens as f64;
2271 return known && fits && observed && pool.latency_n as usize >= LATENCY_MIN;
2272 }
2273 tokens == baseline || known || lambda > 0.0
2274}
2275
2276struct PoolDraws {
2279 grid: Vec<usize>,
2280 shared: Vec<Vec<f64>>,
2281}
2282
2283impl PoolDraws {
2284 fn of(pool: &PoolFit) -> PoolDraws {
2285 let grid = pool_grid();
2286 let shared = grid.iter().map(|&t| (0..DRAWS).map(|k| pool_value(pool, t, mix(hash(&[b"pool-gate"]), k as u64))).collect()).collect();
2287 PoolDraws { grid, shared }
2288 }
2289
2290 fn choose(&self, play: &[f64], feasible: &[bool], lambda: f64) -> usize {
2295 let baseline = crate::whole::TOKENS;
2296 let Some(b) = self.grid.iter().position(|&t| t == baseline) else { return baseline };
2297 let admits = |i: usize| -> bool {
2298 let cost = lambda * (self.grid[i] as f64 - baseline as f64);
2299 let mut gain: Vec<f64> = self.shared[i].iter().zip(&self.shared[b]).map(|(x, y)| x - y - cost).collect();
2300 let q = lower_quantile(&mut gain);
2301 i == b || if self.grid[i] > baseline { q > 0.0 } else { q >= 0.0 }
2302 };
2303 (0..self.grid.len())
2304 .filter(|&i| feasible[i] && admits(i))
2305 .map(|i| (play[i] - lambda * self.grid[i] as f64, self.grid[i]))
2306 .fold(None, |best: Option<(f64, usize)>, x| if best.is_none_or(|b| x.0 > b.0) { Some(x) } else { best })
2307 .map_or(baseline, |(_, t)| t)
2308 }
2309
2310 fn play(&self, pool: &PoolFit, seed: u64) -> Vec<f64> {
2312 self.grid.iter().map(|&t| pool_value(pool, t, seed)).collect()
2313 }
2314}
2315
2316#[must_use]
2323pub fn pool_tokens(fit: &Fit, active: bool, limit_ms: Option<u64>, seed: u64) -> usize {
2324 let baseline = crate::whole::TOKENS;
2325 if !active || fit.pool.yields.is_empty() && fit.lambda == 0.0 {
2326 return baseline;
2327 }
2328 let draws = PoolDraws::of(&fit.pool);
2329 let feasible: Vec<bool> = draws.grid.iter().map(|&t| pool_feasible(&fit.pool, t, limit_ms, fit.lambda)).collect();
2330 draws.choose(&draws.play(&fit.pool, seed), &feasible, fit.lambda)
2331}
2332
2333#[must_use]
2342pub fn price(s: &Summary, fit: &Fit, budget: f64) -> f64 {
2343 let sessions: Vec<(&str, u32)> = s.sessions.iter().rev().take(FIT_SESSIONS).filter(|x| x.checks > 0).map(|x| (x.id.as_str(), x.checks)).collect();
2344 if sessions.is_empty() || budget <= 0.0 {
2345 return 0.0;
2346 }
2347 let draws = PoolDraws::of(&fit.pool);
2348 let plays: Vec<Vec<f64>> = sessions.iter().map(|(id, _)| draws.play(&fit.pool, session_seed(id))).collect();
2349 let spend = |lambda: f64| -> f64 {
2350 let feasible: Vec<bool> = draws.grid.iter().map(|&t| pool_feasible(&fit.pool, t, Some(u64::MAX), lambda)).collect();
2352 sessions.iter().zip(&plays).map(|((_, n), play)| f64::from(*n) * draws.choose(play, &feasible, lambda) as f64).sum::<f64>() / sessions.len() as f64
2353 };
2354 if spend(0.0) <= budget {
2355 return 0.0;
2356 }
2357 let mut hi = 1.0 / crate::whole::TOKENS as f64;
2359 for _ in 0..PRICE_STEPS {
2360 if spend(hi) <= budget {
2361 break;
2362 }
2363 hi *= 2.0;
2364 }
2365 let mut lo = 0.0;
2366 for _ in 0..PRICE_STEPS {
2367 let mid = 0.5 * (lo + hi);
2368 if spend(mid) <= budget { hi = mid } else { lo = mid }
2369 }
2370 hi
2371}
2372
2373#[cfg(test)]
2374mod tests {
2375 use super::*;
2376 use crate::outcome::{self, Outcome};
2377 use crate::select::summary::{self, Summary};
2378
2379 fn unit(lesson: &str, grade: u8, rank: usize) -> Unit<'_> {
2380 Unit { lesson, grade, rank }
2381 }
2382
2383 fn tables(n: usize) -> Vec<Vec<u8>> {
2385 let mut out = vec![Vec::new()];
2386 for len in 1..=n {
2387 out.extend((0..3usize.pow(len as u32)).map(|mut x| {
2388 (0..len)
2389 .map(|_| {
2390 let g = (x % 3) as u8 + 1;
2391 x /= 3;
2392 g
2393 })
2394 .collect()
2395 }));
2396 }
2397 out
2398 }
2399
2400 #[derive(Clone, Copy)]
2403 struct Truth {
2404 b: [f64; P],
2405 lift: f64,
2406 }
2407
2408 impl Truth {
2409 fn p(&self, f: &Feat, route: Route, m: usize) -> f64 {
2410 sigmoid(dot(&self.b, &design(f, route, m)) + if f.stop && f.three { self.lift } else { 0.0 })
2411 }
2412 }
2413
2414 struct Sim {
2418 rng: Rng,
2419 lines: Vec<Outcome>,
2420 ms: i64,
2421 }
2422
2423 type Router<'a> = dyn FnMut(&[u8], bool, &mut Rng) -> (Vec<Route>, Vec<f64>, Vec<f64>, &'static str) + 'a;
2424
2425 impl Sim {
2426 fn new(seed: u64) -> Sim {
2427 Sim { rng: Rng::new(seed), lines: Vec::new(), ms: 0 }
2428 }
2429
2430 fn session(&mut self, session: &str, checks: usize, lessons: usize, truth: &Truth, route: &mut Router) {
2431 for c in 0..checks {
2432 let n = 1 + (self.rng.next_u64() % 4) as usize;
2433 let stop = c % 4 == 3;
2434 let grades: Vec<u8> = (0..n).map(|_| 1 + (self.rng.next_u64() % 3) as u8).collect();
2435 let ids: Vec<String> = (0..n).map(|_| format!("L{}", self.rng.next_u64() % lessons as u64)).collect();
2436 self.check(session, stop, &grades, &ids, truth, route);
2437 }
2438 self.ms += 1;
2439 let asked: Vec<String> =
2440 self.lines.iter().filter(|o| o.session == session && o.kind == outcome::DELIVERY).flat_map(|o| o.lessons.clone()).collect();
2441 if let Some(l) = asked.first() {
2442 self.lines.push(Outcome {
2443 kind: outcome::NOT_CORRECTED.into(),
2444 session: session.into(),
2445 lessons: vec![l.clone()],
2446 ms: self.ms,
2447 reason: outcome::UNJUDGED_NONE.into(),
2448 ..Default::default()
2449 });
2450 }
2451 }
2452
2453 fn check(&mut self, session: &str, stop: bool, grades: &[u8], ids: &[String], truth: &Truth, route: &mut Router) {
2454 self.ms += 1;
2455 let check = crate::select::check_id(session, if stop { "stop" } else { "live" }, self.ms, 1);
2456 let (routes, p, pf, policy) = route(grades, stop, &mut self.rng);
2457 let pick = |r: Route| ids.iter().zip(&routes).filter(|(_, x)| **x == r).map(|(l, _)| l.clone()).collect::<Vec<_>>();
2458 self.lines.push(Outcome {
2459 kind: outcome::GRADED.into(),
2460 session: session.into(),
2461 lessons: ids.to_vec(),
2462 ms: self.ms,
2463 channel: if stop { "stop" } else { "live" }.into(),
2464 grades: grades.to_vec(),
2465 ranks: (0..ids.len() as u16).collect(),
2466 applied: ids.len() as u16,
2467 chosen: pick(Route::Ask),
2468 free: pick(Route::Free),
2469 placebo: pick(Route::Placebo),
2470 psi: 1.0,
2471 p,
2472 pf,
2473 check: check.clone(),
2474 policy: policy.into(),
2475 ..Default::default()
2476 });
2477 let three = grades.contains(&3);
2478 let m = |r: Route| routes.iter().filter(|x| **x == r).count();
2479 for (j, (l, r)) in ids.iter().zip(&routes).enumerate() {
2480 if *r == Route::Not {
2481 continue;
2482 }
2483 self.ms += 1;
2484 let delivered = match r {
2485 Route::Ask => outcome::ROUTE_ASKED,
2486 Route::Free => outcome::ROUTE_FREE,
2487 _ => outcome::ROUTE_PLACEBO,
2488 };
2489 self.lines.push(Outcome {
2490 kind: outcome::DELIVERY.into(),
2491 session: session.into(),
2492 lessons: vec![l.clone()],
2493 ms: self.ms,
2494 check: check.clone(),
2495 route: delivered.into(),
2496 ..Default::default()
2497 });
2498 let f = Feat { grade: grades[j], rank: j, stop, three, pool: 0.0 };
2499 let as_route = if *r == Route::Placebo { Route::Not } else { *r };
2500 let y = self.rng.unit() < truth.p(&f, as_route, m(*r).max(1));
2501 let channel = match r {
2502 Route::Ask => {
2503 if stop {
2504 "stop"
2505 } else {
2506 "live"
2507 }
2508 }
2509 Route::Free => "queued",
2510 _ => outcome::PLACEBO,
2511 };
2512 self.lines.push(Outcome {
2513 kind: if y { outcome::CHANGED } else { outcome::IGNORED }.into(),
2514 session: session.into(),
2515 lessons: vec![l.clone()],
2516 ms: self.ms,
2517 channel: channel.into(),
2518 check: check.clone(),
2519 judge: 2,
2520 ..Default::default()
2521 });
2522 }
2523 }
2524
2525 fn text(&self) -> String {
2526 self.lines.iter().map(|o| format!("{}\n", serde_json::to_string(o).unwrap())).collect()
2527 }
2528
2529 fn summary(&self) -> Summary {
2530 summary::fold(&self.text(), Summary::default())
2531 }
2532 }
2533
2534 fn uniform(grades: &[u8], stop: bool, rng: &mut Rng) -> (Vec<Route>, Vec<f64>, Vec<f64>, &'static str) {
2536 let routes: Vec<Route> = grades
2537 .iter()
2538 .map(|g| if stop && *g == 3 { Route::Ask } else { [Route::Ask, Route::Free, Route::Placebo][(rng.next_u64() % 3) as usize] })
2539 .collect();
2540 let p = grades.iter().map(|g| if stop && *g == 3 { 1.0 } else { 0.3333 }).collect();
2541 let pf = grades.iter().map(|g| if stop && *g == 3 { 0.0 } else { 0.3333 }).collect();
2542 (routes, p, pf, "fixed")
2543 }
2544
2545 fn fixed_rule(grades: &[u8], _: bool, _: &mut Rng) -> (Vec<Route>, Vec<f64>, Vec<f64>, &'static str) {
2547 let ids: Vec<String> = (0..grades.len()).map(|i| format!("x{i}")).collect();
2548 let units: Vec<Unit> = ids.iter().zip(grades).enumerate().map(|(r, (l, g))| unit(l, *g, r)).collect();
2549 let routes = fixed_routes(&units, &[]);
2550 let p = routes.iter().map(|r| f64::from(u8::from(*r == Route::Ask))).collect();
2551 let pf = routes.iter().map(|r| f64::from(u8::from(*r == Route::Free))).collect();
2552 (routes, p, pf, "fixed")
2553 }
2554
2555 fn truth() -> Truth {
2556 let mut b = [0.0; P];
2557 b[B0] = -1.0;
2558 b[G3] = 1.2;
2559 b[ASK] = 0.6;
2560 b[NOT] = -0.8;
2561 Truth { b, lift: 0.0 }
2562 }
2563
2564 fn planted(sessions: usize, checks: usize, seed: u64) -> Summary {
2565 let mut sim = Sim::new(seed);
2566 for i in 0..sessions {
2567 sim.session(&format!("s{i}"), checks, 30, &truth(), &mut uniform);
2568 }
2569 sim.summary()
2570 }
2571
2572 #[test]
2573 fn the_fit_recovers_a_planted_world() {
2574 let s = planted(24, 40, 1);
2575 let fit = fit(&s, None);
2576 let m = fit.model.as_ref().expect("a model");
2577 assert!(m.converged, "{m:?}");
2578 assert!(fit.labels > 600, "{}", fit.labels);
2579 let cov = |k: usize| (0..=k).map(|j| m.chol[k * P + j].powi(2)).sum::<f64>().sqrt();
2580 for k in [B0, G3, ASK, NOT] {
2581 let (got, want, sd) = (m.beta[k], truth().b[k], cov(k));
2582 assert!((got - want).abs() < 3.0 * sd + 0.15, "{}: {got:.3} vs {want} (sd {sd:.3})", NAMES[k]);
2583 }
2584 assert!(m.beta[ASK] > 0.0 && m.beta[NOT] < 0.0, "{:?}", m.beta);
2585 assert!(m.sd[RANK] < ARD_START, "{}", m.sd[RANK]);
2587 assert!(m.lessons.len() <= 30 && m.sessions.len() == 24);
2588 assert!(m.beta.iter().chain(&m.chol).all(|x| x.is_finite()));
2589 }
2590
2591 #[test]
2595 fn a_dated_prior_keeps_its_sd_where_the_labels_say_nothing() {
2596 let mut sim = Sim::new(12);
2597 for i in 0..16 {
2598 sim.session(&format!("s{i}"), 16, 12, &truth(), &mut fixed_rule);
2599 }
2600 let f = fit(&sim.summary(), None);
2601 let m = f.model.as_ref().expect("a model");
2602 for k in [ASK_M, FREE_M] {
2603 assert_eq!(m.sd[k], DILUTION_PRIOR.1, "{}: {:?}", NAMES[k], m.sd);
2604 }
2605 assert_eq!((m.sd[ASK], m.sd[NOT]), (ASK_PRIOR.1, NOT_PRIOR.1));
2606 assert!([RANK, ASK_STOP, ASK_G3, NOT_G, POOL].iter().any(|&k| m.sd[k] != ARD_START), "EM ran on the rest: {:?}", m.sd);
2607 }
2608
2609 #[test]
2613 fn the_pruned_prior_is_the_pruned_models_before_any_label() {
2614 let (full, pruned) = (prior(), prior_pruned());
2615 let (fm, pm) = (full.model.as_ref().expect("a prior"), pruned.model.as_ref().expect("a pruned prior"));
2616 let spec = coefs();
2617 for (k, name) in NAMES.iter().enumerate() {
2618 let want = if PRUNED_KEEPS.contains(&k) { (spec[k].mean, spec[k].sd) } else { (0.0, SD_FLOOR) };
2619 assert!((pm.beta[k] - want.0).abs() < 1e-9 && pm.sd[k] == want.1, "{name}: {} sd {}", pm.beta[k], pm.sd[k]);
2620 }
2621 assert_eq!((pm.spread, fm.spread), ([0.0; FAMILIES], [LEVEL_START; FAMILIES]));
2622 let f = Feat { grade: 2, rank: 3, stop: true, three: false, pool: 0.5 };
2623 let at = |fit: &Fit| predict(fit, &f, Route::Not, 1, "unseen", "L99").expect("a prediction");
2625 assert_ne!(at(&pruned), at(&full), "the full prior's rank, pool and level variance are not the pruned model's");
2626 assert_eq!(pruned.rates, full.rates);
2627 }
2628
2629 #[test]
2634 fn the_pruned_model_keeps_grade_route_and_dilution_alone() {
2635 let mut t = truth();
2636 t.b[RANK] = -0.6;
2637 t.b[STOP] = 0.8;
2638 let mut sim = Sim::new(4);
2639 for i in 0..20 {
2640 sim.session(&format!("s{i}"), 30, 20, &t, &mut uniform);
2641 }
2642 let s = sim.summary();
2643 let (full, pruned) = (fit(&s, None), fit_pruned(&s, None));
2644 let (fm, pm) = (full.model.as_ref().expect("a full model"), pruned.model.as_ref().expect("a pruned model"));
2645 assert!(fm.beta[RANK] < -0.2 && fm.beta[STOP] > 0.2 && !fm.sessions.is_empty(), "the scenario: {:?}", fm.beta);
2646 assert!(pm.current() && pm.converged && pruned.labels == full.labels && pruned.key == full.key);
2647 for (k, name) in NAMES.iter().enumerate() {
2648 if PRUNED_KEEPS.contains(&k) {
2649 assert!(pm.sd[k] > SD_FLOOR, "{name} is kept");
2650 } else {
2651 assert_eq!(pm.sd[k], SD_FLOOR, "{name}");
2652 assert!(pm.beta[k].abs() < 3.0 * SD_FLOOR, "{name}: {}", pm.beta[k]);
2653 }
2654 }
2655 assert!(pm.beta[ASK] > 0.0 && pm.beta[NOT] < 0.0, "{:?}", pm.beta);
2656 assert!(pm.sessions.is_empty() && pm.lessons.is_empty() && pm.judges.is_empty() && pm.judge.is_empty());
2657 assert_eq!(pm.spread, [0.0; FAMILIES]);
2658 assert!(HALF_LIVES.contains(&pruned.half_life));
2659 let f = Feat { grade: 2, rank: 3, stop: true, three: false, pool: 0.0 };
2660 let at = |session: &str, lesson: &str| predict(&pruned, &f, Route::Free, 2, session, lesson);
2661 assert_eq!(at("s0", "L1"), at("unseen", "L99"), "no level effect");
2662 assert_ne!(predict(&full, &f, Route::Free, 2, "s0", "L1"), predict(&full, &f, Route::Free, 2, "unseen", "L99"), "the full model's levels");
2663 }
2664
2665 #[test]
2669 fn eliminating_a_family_equals_solving_it_densely() {
2670 let mut rng = Rng::new(9);
2671 let rows: Vec<([f64; 3], usize, f64)> = (0..400)
2672 .map(|i| {
2673 let x = [1.0, rng.unit() * 2.0 - 1.0, f64::from(u8::from(rng.unit() < 0.4))];
2674 let l = i % 4;
2675 let eta = -0.3 + 0.8 * x[1] - 0.5 * x[2] + [0.6, -0.4, 0.1, 0.0][l];
2676 (x, l, f64::from(u8::from(rng.unit() < sigmoid(eta))))
2677 })
2678 .collect();
2679 let spread = 0.7;
2680 let spec3 = [Coef { mean: 0.0, sd: 2.5, ard: false, nonpos: false }; 3];
2681 let obs: Vec<Obs<3, 1>> = rows.iter().map(|(x, l, y)| Obs { x: *x, y: *y, w: 1.0 + f64::from(*l as u8) * 0.1, lv: [Some(*l)] }).collect();
2682 let a = map(&obs, [4], &spec3, &[2.5; 3], &[spread], None).unwrap();
2683 let mut spec7 = [Coef { mean: 0.0, sd: spread, ard: false, nonpos: false }; 7];
2684 spec7[..3].copy_from_slice(&spec3);
2685 let sd7: [f64; 7] = std::array::from_fn(|k| spec7[k].sd);
2686 let dense: Vec<Obs<7, 0>> = rows
2687 .iter()
2688 .map(|(x, l, y)| {
2689 let mut z = [0.0; 7];
2690 z[..3].copy_from_slice(x);
2691 z[3 + l] = 1.0;
2692 Obs { x: z, y: *y, w: 1.0 + f64::from(*l as u8) * 0.1, lv: [] }
2693 })
2694 .collect();
2695 let b = map(&dense, [], &spec7, &sd7, &[], None).unwrap();
2696 assert!(a.converged && b.converged);
2697 for k in 0..3 {
2698 assert!((a.beta[k] - b.beta[k]).abs() < 1e-6, "{k}: {} {}", a.beta[k], b.beta[k]);
2699 for j in 0..3 {
2700 assert!((a.cov[k][j] - b.cov[k][j]).abs() < 1e-6, "cov {k} {j}");
2701 }
2702 }
2703 for l in 0..4 {
2704 assert!((a.u[0][l] - b.beta[3 + l]).abs() < 1e-6, "level {l}");
2705 assert!((a.var_u[0][l] - b.cov[3 + l][3 + l]).abs() < 1e-6, "var {l}: {} {}", a.var_u[0][l], b.cov[3 + l][3 + l]);
2706 }
2707 }
2708
2709 #[test]
2710 fn a_bounded_slope_stops_at_zero_and_a_negative_one_is_learned() {
2711 let mut rng = Rng::new(3);
2712 let obs: Vec<Obs<3, 0>> = (0..600)
2713 .map(|_| {
2714 let x = [1.0, rng.unit() * 3.0, rng.unit() * 3.0];
2715 let eta = -0.5 + 0.9 * x[1] - 0.7 * x[2];
2716 Obs { x, y: f64::from(u8::from(rng.unit() < sigmoid(eta))), w: 1.0, lv: [] }
2717 })
2718 .collect();
2719 let bounded = Coef { mean: 0.0, sd: 2.5, ard: false, nonpos: true };
2720 let spec = [Coef { nonpos: false, ..bounded }, bounded, bounded];
2721 let e = map(&obs, [], &spec, &[2.5; 3], &[], None).unwrap();
2722 assert_eq!(e.beta[1], 0.0, "a positive effect under a ≤ 0 bound sits at the bound");
2723 assert!(e.beta[2] < -0.3, "{:?}", e.beta);
2724 }
2725
2726 fn brute(v: &Valuer, nofree: &[bool], rail: usize, cap: usize) -> f64 {
2728 let n = v.n();
2729 let mut best = f64::NEG_INFINITY;
2730 for code in 0..3usize.pow(n as u32) {
2731 let mut x = code;
2732 let routes: Vec<Route> = (0..n)
2733 .map(|_| {
2734 let r = [Route::Ask, Route::Free, Route::Not][x % 3];
2735 x /= 3;
2736 r
2737 })
2738 .collect();
2739 let (m_q, m_f) = (routes.iter().filter(|r| **r == Route::Ask).count(), routes.iter().filter(|r| **r == Route::Free).count());
2740 let forced_ok = v.feats.iter().zip(&routes).all(|(f, r)| !forced(f) || *r == Route::Ask);
2741 let free_ok = nofree.iter().zip(&routes).all(|(x, r)| !x || *r != Route::Free);
2742 let f = v.feats.iter().filter(|f| forced(f)).count();
2743 if forced_ok && free_ok && m_q <= rail.max(f) && m_f <= cap {
2744 best = best.max(v.value(&routes));
2745 }
2746 }
2747 best
2748 }
2749
2750 #[test]
2751 fn the_solver_equals_brute_force() {
2752 let mut rng = Rng::new(11);
2753 for case in 0..400 {
2754 let n = case % 8;
2755 let stop = case % 3 == 0;
2756 let feats: Vec<Feat> = (0..n).map(|r| Feat { grade: 1 + (rng.next_u64() % 3) as u8, rank: r, stop, three: false, pool: 0.0 }).collect();
2757 let mut b = [0.0; P];
2758 for x in &mut b {
2759 *x = 2.0 * normal(&mut rng);
2760 }
2761 b[ASK_M] = -rng.unit() * 2.0;
2762 b[FREE_M] = -rng.unit() * 2.0;
2763 let t = Theta { b, base: (0..n).map(|_| normal(&mut rng)).collect() };
2764 let rates = Rates { free_stop: rng.unit(), free_live: rng.unit(), later: rng.unit(), wait: 0.0, ask_live: 1.0 };
2765 let v = Valuer { t: &t, feats: &feats, rates: &rates };
2766 let nofree: Vec<bool> = (0..n).map(|_| case % 2 == 1 && rng.unit() < 0.4).collect();
2768 let table = Table::of(&v, &nofree);
2769 for (rail, cap) in [(20, 10), (2, 1), (1, 3), (0, 0)] {
2770 let got = solve(&table, rail, cap);
2771 let want = brute(&v, &nofree, rail, cap);
2772 assert!((v.value(&got) - want).abs() < 1e-12, "case {case} rail {rail} cap {cap}: {} vs {want} {got:?}", v.value(&got));
2773 assert!((table.value(&got) - v.value(&got)).abs() < 1e-12, "the table values as the model does");
2774 assert!(feats.iter().zip(&got).all(|(f, r)| !forced(f) || *r == Route::Ask));
2775 assert!(nofree.iter().zip(&got).all(|(x, r)| !x || *r != Route::Free), "case {case}: {nofree:?} {got:?}");
2776 }
2777 }
2778 }
2779
2780 #[test]
2781 fn the_value_of_a_route_is_the_models_lift() {
2782 let f = Feat { grade: 3, rank: 2, stop: false, three: true, pool: 0.3 };
2783 let mut b = [0.0; P];
2784 for (i, x) in b.iter_mut().enumerate() {
2785 *x = 0.1 * i as f64 - 0.6;
2786 }
2787 let base = dot(&b, &design(&f, Route::Free, 1)) + 0.25;
2788 let t = Theta { b, base: vec![base] };
2789 let rates = Rates { free_stop: 0.5, free_live: 0.8, later: 0.4, wait: 0.0, ask_live: 0.6 };
2790 let v = Valuer { t: &t, feats: std::slice::from_ref(&f), rates: &rates };
2791 let p = |r: Route, m: usize| sigmoid(dot(&b, &design(&f, r, m)) + 0.25);
2792 assert!((v.ask(0, 3) - 0.6 * (p(Route::Ask, 3) - p(Route::Not, 1))).abs() < 1e-12);
2794 assert!((v.free(0, 2) - 0.8 * (p(Route::Free, 2) - p(Route::Not, 1))).abs() < 1e-12);
2795 for r in [Route::Not, Route::Placebo] {
2798 assert_eq!(v.value(&[r]), 0.0, "during the turn a lesson left for later counts nothing here: {r:?}");
2799 }
2800 let stop = Feat { stop: true, ..f };
2801 let t = Theta { b, base: vec![dot(&b, &design(&stop, Route::Free, 1)) + 0.25] };
2802 let v = Valuer { t: &t, feats: std::slice::from_ref(&stop), rates: &rates };
2803 let p = |r: Route, m: usize| sigmoid(dot(&b, &design(&stop, r, m)) + 0.25);
2804 assert!((v.ask(0, 3) - (p(Route::Ask, 3) - p(Route::Not, 1))).abs() < 1e-12, "a stop's question blocks: it arrives");
2805 assert!((v.free(0, 2) - 0.5 * (p(Route::Free, 2) - p(Route::Not, 1))).abs() < 1e-12);
2806 for r in [Route::Not, Route::Placebo] {
2807 assert_eq!(v.value(&[r]), 0.0, "at the stop a lesson not delivered counts nothing: {r:?}");
2808 }
2809 }
2810
2811 #[test]
2812 fn a_live_question_without_a_3_arrives_as_a_hand_over_does() {
2813 let mut b = [0.0; P];
2814 for (i, x) in b.iter_mut().enumerate() {
2815 *x = 0.1 * i as f64 - 0.6;
2816 }
2817 let t = Theta { b, base: vec![0.25] };
2818 let rates = Rates { free_stop: 0.5, free_live: 0.3, later: 0.4, wait: 0.0, ask_live: 0.9 };
2819 let ask = |three: bool| {
2820 let f = Feat { grade: 2, rank: 0, stop: false, three, pool: 0.0 };
2821 Valuer { t: &t, feats: std::slice::from_ref(&f), rates: &rates }.ask(0, 1)
2822 };
2823 let lift = sigmoid(0.25 + b[ASK]) - sigmoid(0.25 + b[NOT]);
2825 assert!(lift.abs() > 0.01);
2826 assert!((ask(true) - 0.9 * lift).abs() < 1e-12);
2828 assert!((ask(false) - 0.3 * lift).abs() < 1e-12);
2831 }
2832
2833 fn open_gate() -> Gate {
2835 Gate { gain: vec![0.0; DRAWS], allowance: 5.0, ..Gate::default() }
2836 }
2837
2838 #[test]
2839 fn cold_start_is_the_fixed_rule() {
2840 let fit_model = fit(&planted(6, 20, 2), None);
2841 assert!(fit_model.model.is_some());
2842 let wrong_len = Gate { gain: vec![0.0; 3], allowance: 1.0, ..Gate::default() };
2843 let ids: Vec<String> = (0..14).map(|i| format!("L{i}")).collect();
2844 let mut grades = tables(5);
2845 grades.extend([vec![3; 12], [vec![1; 11], vec![2; 3]].concat()]);
2846 for g in &grades {
2847 let units: Vec<Unit> = ids.iter().zip(g).enumerate().map(|(r, (l, g))| unit(l, *g, r)).collect();
2848 let want = fixed_routes(&units, &[]);
2849 for (fit, gate) in [(&Fit::default(), &open_gate()), (&fit_model, &Gate::default()), (&fit_model, &wrong_len)] {
2850 for stop in [true, false] {
2851 let d = choose(&units, &Check { session: "s", stop, tokens: 0 }, fit, gate, ALPHA, 7, &[]);
2852 assert_eq!(d.routes, want, "{g:?}");
2853 assert_eq!(d.p, want.iter().map(|r| f64::from(u8::from(*r == Route::Ask))).collect::<Vec<_>>());
2854 assert_eq!(d.pf, want.iter().map(|r| f64::from(u8::from(*r == Route::Free))).collect::<Vec<_>>());
2855 assert_eq!((d.source, d.gate), ("fixed", 0.0));
2856 }
2857 }
2858 }
2859 }
2860
2861 #[test]
2862 fn n_zero_and_n_one() {
2863 let empty = fit(&Summary::default(), None);
2864 assert!(empty.model.is_none() && empty.labels == 0);
2865 let mut sim = Sim::new(5);
2866 let mut one = |_: &[u8], _: bool, _: &mut Rng| (vec![Route::Ask], vec![1.0], vec![0.0], "fixed");
2867 sim.check("s", false, &[2], &["A".to_string()], &truth(), &mut one);
2868 let f = fit(&sim.summary(), None);
2869 let m = f.model.as_ref().expect("one label is a model");
2870 assert_eq!(f.labels, 1);
2871 assert!(m.beta.iter().chain(&m.chol).all(|x| x.is_finite()), "{m:?}");
2872 let gate = open_gate();
2873 assert!(choose(&[], &Check { session: "s", stop: true, tokens: 0 }, &f, &gate, ALPHA, 1, &[]).routes.is_empty());
2874 for g in 1..=3 {
2875 let d = choose(&[unit("A", g, 0)], &Check { session: "s", stop: true, tokens: 0 }, &f, &gate, ALPHA, 1, &[]);
2876 assert_eq!(d.routes.len(), 1);
2877 assert!((d.p[0] + d.pf[0]) <= 1.0 + 1e-9);
2878 }
2879 }
2880
2881 #[test]
2884 fn decisions_are_deterministic_and_only_the_test_varies() {
2885 let f = fit(&planted(8, 20, 3), None);
2886 let gate = open_gate();
2887 let units = [unit("L1", 2, 0), unit("L2", 1, 1), unit("L3", 3, 2)];
2888 let check = Check { session: "new", stop: false, tokens: 0 };
2889 let a = choose(&units, &check, &f, &gate, ALPHA, 42, &[]);
2890 assert_eq!(a, choose(&units, &check, &f, &gate, ALPHA, 42, &[]));
2891 let plays: std::collections::HashSet<Vec<Route>> = (0..60).map(|s| choose(&units, &check, &f, &gate, ALPHA, s, &[]).routes).collect();
2892 assert!(plays.len() > 1, "the test's count varies with the check: {plays:?}");
2893 let stop = Check { stop: true, ..check };
2894 let at_stop: std::collections::HashSet<Vec<Route>> = (0..60).map(|s| choose(&units, &stop, &f, &gate, ALPHA, s, &[]).routes).collect();
2895 assert_eq!(at_stop.len(), 1, "no Thompson play: {at_stop:?}");
2896 }
2897
2898 #[test]
2899 fn propensities_match_the_frequency_of_the_play() {
2900 let f = fit(&planted(6, 15, 4), None);
2901 let gate = open_gate();
2902 let units = [unit("L1", 2, 0), unit("L2", 1, 1), unit("L3", 1, 2)];
2903 let check = Check { session: "new", stop: false, tokens: 0 };
2904 let n = 300;
2905 let (mut asked, mut free) = ([0.0; 3], [0.0; 3]);
2906 let (mut p, mut pf) = ([0.0; 3], [0.0; 3]);
2907 for seed in 0..n {
2908 let d = choose(&units, &check, &f, &gate, ALPHA, seed, &[]);
2909 for j in 0..3 {
2910 asked[j] += f64::from(u8::from(d.routes[j] == Route::Ask)) / n as f64;
2911 free[j] += f64::from(u8::from(d.routes[j] == Route::Free)) / n as f64;
2912 p[j] += d.p[j] / n as f64;
2913 pf[j] += d.pf[j] / n as f64;
2914 }
2915 }
2916 for j in 0..3 {
2917 for (freq, prop) in [(asked[j], p[j]), (free[j], pf[j])] {
2918 let se = (prop * (1.0 - prop) / n as f64).sqrt();
2919 assert!((freq - prop).abs() <= 4.0 * se + 0.01, "unit {j}: played {freq:.3}, propensity {prop:.3}");
2920 }
2921 }
2922 assert!(asked.iter().chain(&free).any(|x| *x > 0.05 && *x < 0.95), "the play varies: {asked:?} {free:?}");
2923 }
2924
2925 fn free_beats_ask() -> Fit {
2927 let mut sim = Sim::new(8);
2928 let mut t = truth();
2929 t.b[ASK] = -3.0;
2930 for i in 0..10 {
2931 sim.session(&format!("s{i}"), 12, 10, &t, &mut uniform);
2932 }
2933 fit(&sim.summary(), None)
2934 }
2935
2936 #[test]
2939 fn a_stop_grade_three_is_never_routed_away_from_asking() {
2940 let f = free_beats_ask();
2941 let gate = open_gate();
2942 let units = [unit("L1", 3, 0), unit("L2", 1, 1)];
2943 let live = choose(&units, &Check { session: "n", stop: false, tokens: 0 }, &f, &gate, ALPHA, 1, &[]);
2944 assert_eq!(live.routes[0], Route::Free, "during the turn the model hands it over: {live:?}");
2945 for seed in 0..20 {
2946 let d = choose(&units, &Check { session: "n", stop: true, tokens: 0 }, &f, &gate, ALPHA, seed, &[]);
2947 assert_eq!((d.routes[0], d.p[0], d.pf[0]), (Route::Ask, 1.0, 0.0), "{d:?}");
2948 }
2949 }
2950
2951 #[test]
2955 fn a_lesson_the_session_had_never_goes_free() {
2956 let f = free_beats_ask();
2957 let gate = open_gate();
2958 let units = [unit("L1", 2, 0), unit("L2", 1, 1), unit("L3", 1, 2)];
2959 let had = ["L2".to_string()];
2960 let check = Check { session: "n", stop: false, tokens: 0 };
2961 assert_eq!(fixed_routes(&units, &had), [Route::Free, Route::Not, Route::Free]);
2962 let cold = choose(&units, &check, &Fit::default(), &gate, ALPHA, 1, &had);
2963 assert_eq!((cold.routes, cold.pf[1], cold.source), (fixed_routes(&units, &had), 0.0, "fixed"));
2964 let mut freed = 0;
2965 for seed in 0..20 {
2966 freed += usize::from(choose(&units, &check, &f, &gate, ALPHA, seed, &[]).routes[1] == Route::Free);
2967 let d = choose(&units, &check, &f, &gate, ALPHA, seed, &had);
2968 assert!(d.routes[1] != Route::Free && d.pf[1] == 0.0, "{d:?}");
2969 }
2970 assert!(freed > 0, "the scenario: without `had` the model hands L2 over");
2971 let block = [unit("L2", 3, 0), unit("L1", 1, 1)];
2973 let d = choose(&block, &Check { stop: true, ..check }, &f, &gate, ALPHA, 3, &had);
2974 assert_eq!((d.routes[0], d.p[0]), (Route::Ask, 1.0));
2975 }
2976
2977 #[test]
2980 fn a_check_past_the_solvers_size_takes_the_fixed_rule() {
2981 let f = free_beats_ask();
2982 let gate = open_gate();
2983 let ids: Vec<String> = (0..=SOLVE_MAX).map(|i| format!("L{i}")).collect();
2984 let units: Vec<Unit> = ids.iter().enumerate().map(|(r, l)| unit(l, if r == 0 { 3 } else { 2 }, r)).collect();
2986 let check = Check { session: "n", stop: false, tokens: 0 };
2987 let at = choose(&units[..SOLVE_MAX], &check, &f, &gate, ALPHA, 1, &[]);
2988 assert_eq!(at.source, EVIDENCE, "{at:?}");
2989 assert_ne!(at.routes, fixed_routes(&units[..SOLVE_MAX], &[]), "at the size the model hands them over: {at:?}");
2990 let past = choose(&units, &check, &f, &gate, ALPHA, 1, &[]);
2991 let base = fixed_routes(&units, &[]);
2992 assert_eq!((&past.routes, past.source, past.risk), (&base, "size", 0.0));
2993 assert_eq!(past.p, base.iter().map(|r| f64::from(u8::from(*r == Route::Ask))).collect::<Vec<_>>(), "its choice is its propensity");
2994 }
2995
2996 #[test]
2999 fn a_confounded_history_keeps_the_stop_three_asked() {
3000 let mut sim = Sim::new(12);
3001 for i in 0..12 {
3002 sim.session(&format!("s{i}"), 16, 12, &truth(), &mut fixed_rule);
3003 }
3004 let s = sim.summary();
3005 assert!(s.applies.iter().all(|a| a.grade != 3 || a.route == Route::Ask || a.route == Route::Not), "never free");
3006 let f = fit(&s, None);
3007 let gate = open_gate();
3008 for seed in 0..10 {
3009 let d = choose(&[unit("L1", 3, 0), unit("L2", 2, 1)], &Check { session: "n", stop: true, tokens: 0 }, &f, &gate, ALPHA, seed, &[]);
3010 assert_eq!((d.routes[0], d.p[0]), (Route::Ask, 1.0));
3011 }
3012 }
3013
3014 fn three_stop_world(lift: f64) -> Summary {
3017 let mut t = truth();
3018 t.b[ASK] = 0.0;
3019 t.b[NOT] = 0.0;
3020 t.lift = lift;
3021 let mut cold = |grades: &[u8], stop: bool, rng: &mut Rng| {
3022 let (routes, p, pf, _) = fixed_rule(grades, stop, rng);
3023 (routes, p, pf, "learned")
3024 };
3025 let mut sim = Sim::new(1);
3026 for i in 0..40 {
3027 sim.session(&format!("s{i}"), 20, 12, &t, &mut cold);
3028 }
3029 summary::with_fit(sim.summary(), None)
3030 }
3031
3032 #[test]
3037 fn a_lift_of_the_stops_with_a_three_is_not_taken_for_asking() {
3038 let (lifted, plain) = (three_stop_world(1.5), three_stop_world(0.0));
3039 let units = [unit("A", 2, 0), unit("B", 1, 1)];
3040 assert_eq!(fixed_routes(&units, &[]), [Route::Free; 2], "the fixed rule hands both over");
3041 let asked = |s: &Summary| -> f64 {
3042 assert!(remaining(&s.gate) > 0.0, "the gate opened on the history: {:?}", s.gate.allowance);
3043 let seeds = 20u32;
3044 (0..seeds)
3045 .map(|seed| {
3046 let d = choose(&units, &Check { session: "new", stop: true, tokens: 0 }, &s.fit, &s.gate, ALPHA, u64::from(seed), &[]);
3047 assert_ne!(d.source, "fixed", "{d:?}");
3048 d.p[0]
3049 })
3050 .sum::<f64>()
3051 / f64::from(seeds)
3052 };
3053 let (with_lift, without) = (asked(&lifted), asked(&plain));
3054 let m = lifted.fit.model.as_ref().expect("a model");
3055 let sd = |k: usize| (0..=k).map(|j| m.chol[k * P + j].powi(2)).sum::<f64>().sqrt();
3056 let (ask_stop, twin) = ((m.beta[ASK_STOP], sd(ASK_STOP)), (m.beta[THREE_STOP], sd(THREE_STOP)));
3057 let seen = format!("asked {with_lift:.3} with the lift, {without:.3} without; ask_stop {ask_stop:.3?}, three_stop {twin:.3?}");
3058 assert!(with_lift < 0.9, "a stop without a 3 is not blocked by the lift: {seen}");
3059 assert!((with_lift - without).abs() < 0.15, "{seen}");
3060 assert!(ask_stop.0.abs() < 2.0 * ask_stop.1, "asking at the stop is not identified: {seen}");
3061 assert!(twin.0 > 0.75, "the lift is the twin's: {seen}");
3062 }
3063
3064 #[test]
3069 fn the_gate_never_risks_more_than_its_alpha_share() {
3070 let spec = coefs();
3073 let mut beta: Vec<f64> = spec.iter().map(|c| c.mean).collect();
3074 beta[ASK] = -0.4;
3075 beta[B0] = -0.5;
3076 let mut chol = vec![0.0; P * P];
3077 for k in 0..P {
3078 chol[k * P + k] = if k == ASK { 1.0 } else { 0.3 };
3079 }
3080 let model = Model { beta, chol, sd: spec.iter().map(|c| c.sd).collect(), spread: [0.3; FAMILIES], ..Model::default() };
3081 let rates = Rates { free_stop: 1.0, free_live: 1.0, later: 0.5, wait: 0.0, ask_live: 1.0 };
3082 let f = Fit { key: "test".into(), model: Some(model), rates, ..Fit::default() };
3083 let ids: Vec<String> = (0..10).map(|i| format!("L{i}")).collect();
3084 let mut refused = 0;
3085 let mut run = |gated: bool, rng: &mut Rng| -> (bool, usize) {
3086 let mut gate = Gate::default();
3087 let mut broke = false;
3088 let mut deviated = 0;
3089 for c in 0..200 {
3090 let n = 1 + (rng.next_u64() % 3) as usize;
3091 let units: Vec<Unit> = (0..n).map(|r| unit(&ids[(rng.next_u64() % 10) as usize], 1 + (rng.next_u64() % 3) as u8, r)).collect();
3092 let session = format!("t{}", c / 10);
3093 let check = Check { session: &session, stop: c % 5 == 0, tokens: 0 };
3094 let routes = if c < 10 {
3095 fixed_routes(&units, &[])
3096 } else if gated {
3097 let d = choose(&units, &check, &f, &gate, ALPHA, c as u64, &[]);
3098 match d.source {
3101 GATE => {
3102 assert_eq!((&d.routes, d.risk), (&fixed_routes(&units, &[]), 0.0), "{d:?}");
3103 refused += 1;
3104 }
3105 _ if d.routes != fixed_routes(&units, &[]) => {
3106 assert_eq!(d.source, EXPERIMENT, "{d:?}");
3107 assert!(d.risk <= d.gate + 1e-12, "{d:?}");
3108 }
3109 _ => {}
3110 }
3111 d.routes
3112 } else {
3113 let open = Gate { allowance: f64::INFINITY, ..gate.clone() };
3114 choose(&units, &check, &f, &open, ALPHA, c as u64, &[]).routes
3115 };
3116 deviated += usize::from(routes != fixed_routes(&units, &[]));
3117 ledger_add_without(&f, &mut gate, &check, &units, &routes, &[]);
3118 broke |= lower_quantile(&mut gate.gain.clone()) < -gate.allowance - 1e-9;
3119 }
3120 (broke, deviated)
3121 };
3122 let (broke, deviated) = run(true, &mut Rng::new(77));
3123 assert!(!broke && deviated > 0, "gated: broke {broke}, deviated {deviated}");
3124 let (broke, _) = run(false, &mut Rng::new(77));
3125 assert!(broke, "ungated, the same run overspends: the scenario needs the gate");
3126 assert!(refused > 0, "the gate refused deviations, logged as `gate`");
3127 }
3128
3129 #[test]
3130 fn no_price_without_a_budget_and_spend_within_one() {
3131 let mut s = planted(10, 12, 6);
3132 s.pool.bins = (0..6).map(|b| (400 * (1u64 << b), 200 / (b as u64 + 1))).collect();
3134 s.pool.shown = 40 * 120;
3135 s.pool.tokens = 12_000 * 120;
3136 let free = fit(&s, None);
3137 assert_eq!((free.lambda, free.budget), (0.0, None), "no budget, no price");
3138 let checks = s.sessions.iter().map(|x| f64::from(x.checks)).sum::<f64>() / s.sessions.len() as f64;
3139 let spend = |f: &Fit| {
3140 s.sessions.iter().map(|x| f64::from(x.checks) * pool_tokens(f, true, None, session_seed(&x.id)) as f64).sum::<f64>() / s.sessions.len() as f64
3141 };
3142 assert!(spend(&free) >= checks * crate::whole::TOKENS as f64 - 1e-6, "unpriced, the pool keeps the baseline");
3143 let budget = checks * 9_000.0;
3144 let priced = fit(&s, Some(budget));
3145 assert!(priced.lambda > 0.0, "a budget below today's spend prices tokens");
3146 assert!(spend(&priced) <= budget + 1e-6, "{} > {budget}", spend(&priced));
3147 let roomy = fit(&s, Some(checks * 1e6));
3148 assert_eq!(roomy.lambda, 0.0, "a budget the baseline meets costs nothing");
3149 }
3150
3151 #[test]
3155 fn a_budget_holds_where_pools_may_grow() {
3156 let mut s = Summary::default();
3157 for i in 0..10 {
3158 s.sessions.push(summary::Session { id: format!("s{i}"), checks: 10, ..summary::Session::default() });
3159 }
3160 let pool = PoolFit {
3164 yields: vec![(50.0, 10.0); 8],
3165 per_lesson: 300.0,
3166 latency: (1000.0, 0.1, 500.0),
3167 latency_n: LATENCY_MIN as u32,
3168 latency_reach: 24_000,
3169 not: (-1.0, 0.1),
3170 ..PoolFit::default()
3171 };
3172 let spend = |f: &Fit, limit: Option<u64>| {
3173 s.sessions.iter().map(|x| f64::from(x.checks) * pool_tokens(f, true, limit, session_seed(&x.id)) as f64).sum::<f64>() / s.sessions.len() as f64
3174 };
3175 let baseline = 10.0 * crate::whole::TOKENS as f64;
3176 let unpriced = Fit { pool: pool.clone(), ..Fit::default() };
3177 assert!(spend(&unpriced, Some(60_000)) > 1.5 * baseline, "the scenario: unpriced, the pools grow");
3178 for budget in [baseline, 1.2 * baseline, 0.9 * baseline] {
3179 let mut f = Fit { pool: pool.clone(), budget: Some(budget), ..Fit::default() };
3180 f.lambda = price(&s, &f, budget);
3181 assert!(f.lambda > 0.0, "budget {budget}: growth is priced");
3182 for limit in [None, Some(3_500), Some(60_000)] {
3183 assert!(spend(&f, limit) <= budget + 1e-6, "budget {budget}, limit {limit:?}: spent {} (λ {:e})", spend(&f, limit), f.lambda);
3184 }
3185 }
3186 }
3187
3188 #[test]
3189 fn the_pool_stays_at_the_baseline_without_evidence_or_a_limit() {
3190 let f = Fit::default();
3191 assert_eq!(pool_tokens(&f, true, Some(60_000), 1), crate::whole::TOKENS);
3192 assert_eq!(pool_tokens(&f, false, Some(60_000), 1), crate::whole::TOKENS);
3193 let mut s = planted(3, 6, 7);
3194 s.pool.bins = (0..6).map(|b| (100 * (1u64 << b), 90)).collect();
3195 s.pool.shown = 40 * 18;
3196 s.pool.tokens = 12_000 * 18;
3197 let f = fit(&s, None);
3198 assert!(pool_tokens(&f, true, Some(60_000), 1) <= crate::whole::TOKENS);
3200 assert_eq!(pool_grid(), [6000, 8485, 12000, 16971, 24000]);
3201 }
3202
3203 #[test]
3207 fn a_larger_pool_is_admitted_only_on_a_positive_gain() {
3208 let grid = pool_grid();
3209 let b = grid.iter().position(|&t| t == crate::whole::TOKENS).unwrap();
3210 let tied = PoolDraws { grid: grid.clone(), shared: vec![vec![1.0; DRAWS]; grid.len()] };
3211 let mut play = vec![0.0; grid.len()];
3212 play[grid.len() - 1] = 5.0;
3213 let feasible: Vec<bool> = (0..grid.len()).map(|i| i >= b).collect();
3215 assert_eq!(tied.choose(&play, &feasible, 0.0), crate::whole::TOKENS, "a tie at 0 keeps the baseline");
3216 let mut ahead = tied.shared.clone();
3217 ahead[grid.len() - 1] = vec![1.5; DRAWS];
3218 let ahead = PoolDraws { grid: grid.clone(), shared: ahead };
3219 assert_eq!(ahead.choose(&play, &feasible, 0.0), grid[grid.len() - 1], "a gain > 0 at every draw admits it");
3220 let mut smaller = vec![0.0; grid.len()];
3221 smaller[0] = 5.0;
3222 assert_eq!(tied.choose(&smaller, &[true; 5], 0.0), grid[0], "a smaller pool on a tie, as before");
3223 }
3224
3225 #[test]
3229 fn a_pool_grows_only_where_its_latency_was_observed() {
3230 let mut s = Summary::default();
3231 s.pool.latency = (0..204u32).map(|i| (if i / 2 < 5 { 9_500 } else { 11_000 + i / 2 }, 16_000 + 3_000 * (i % 2))).collect();
3233 s.pool.bins = vec![(400, 300); 9];
3235 (s.pool.shown, s.pool.tokens) = (100, 9_100);
3236 let mut pool = PoolFit::default();
3237 pool_state(&s, &mut pool);
3238 assert!(pool.yields.len() == 9 && pool.per_lesson == 91.0, "{pool:?}");
3239 assert_eq!((pool.latency.1, pool.latency_n), (0.0, 204), "the scenario: a flat fit, enough calls");
3240 assert!(pool.latency_reach >= 11_000 && pool.latency_reach < 12_000, "{}", pool.latency_reach);
3241 let roomy = Some(u64::MAX);
3242 assert!(pool_feasible(&pool, crate::whole::TOKENS, roomy, 0.0));
3243 for t in pool_grid().into_iter().filter(|&t| t > crate::whole::TOKENS) {
3244 assert!(!pool_feasible(&pool, t, roomy, 0.0), "{t}: no call read near it");
3245 }
3246 s.pool.latency.extend((0..LATENCY_MIN as u32).map(|i| (15_300 + i, 20_000)));
3248 pool_state(&s, &mut pool);
3249 assert_eq!(pool.latency_reach, 15_300);
3250 assert!(pool_feasible(&pool, 16_971, roomy, 0.0) && !pool_feasible(&pool, 24_000, roomy, 0.0));
3251 assert!(!pool_feasible(&PoolFit { latency_reach: 0, ..pool }, 16_971, roomy, 0.0));
3253 }
3254
3255 #[test]
3258 fn the_yields_stop_at_a_bin_no_lesson_was_shown_in() {
3259 let mut s = Summary::default();
3260 s.pool.bins = vec![(40, 20), (80, 30), (0, 0), (3, 1)];
3261 let mut pool = PoolFit::default();
3262 pool_state(&s, &mut pool);
3263 assert_eq!(pool.yields, pava(&[(21.0, 21.0), (31.0, 51.0)]), "the gap ends them");
3264 s.pool.bins[2] = (5, 1);
3265 pool_state(&s, &mut pool);
3266 assert_eq!(pool.yields.len(), 4, "filled, the tail's bin past it counts");
3267 }
3268
3269 #[test]
3270 fn pava_keeps_the_yields_from_rising() {
3271 let y = pava(&[(5.0, 5.0), (8.0, 2.0), (1.0, 9.0), (2.0, 8.0)]);
3272 let means: Vec<f64> = y.iter().map(|(a, b)| a / (a + b)).collect();
3273 assert!(means.windows(2).all(|w| w[0] >= w[1] - 1e-12), "{means:?}");
3274 assert_eq!(y[0], y[1], "the rising pair merged");
3275 assert_eq!(y[0], (13.0, 7.0));
3276 }
3277
3278 #[test]
3279 fn the_samplers_have_their_moments() {
3280 let mut rng = Rng::new(1);
3281 let n = 20_000;
3282 let z: Vec<f64> = (0..n).map(|_| normal(&mut rng)).collect();
3283 let mean = z.iter().sum::<f64>() / n as f64;
3284 let var = z.iter().map(|x| (x - mean).powi(2)).sum::<f64>() / n as f64;
3285 assert!(mean.abs() < 0.03 && (var - 1.0).abs() < 0.05, "{mean} {var}");
3286 for shape in [0.4, 1.0, 3.5] {
3287 let g = (0..n).map(|_| gamma(shape, &mut rng)).sum::<f64>() / n as f64;
3288 assert!((g - shape).abs() < 0.05 * shape.max(1.0), "{shape}: {g}");
3289 }
3290 let b = (0..n).map(|_| beta_draw(2.0, 6.0, &mut rng)).sum::<f64>() / n as f64;
3291 assert!((b - 0.25).abs() < 0.01, "{b}");
3292 }
3293
3294 #[test]
3296 fn the_upper_bound_covers_the_mean_at_every_time() {
3297 let mut rng = Rng::new(21);
3298 let mut missed = 0;
3299 for _ in 0..400 {
3300 let x: Vec<f64> = (0..200).map(|_| normal(&mut rng) * 0.4 - 0.1).collect();
3301 missed += usize::from(upper_bound(&x) < -0.1);
3302 }
3303 assert!(missed <= 400 * 7 / 100, "{missed} of 400 below the mean at some time (δ = 5%)");
3304 assert_eq!(upper_bound(&[-5.0; CS_MIN - 1]), f64::INFINITY, "too few sessions to trip");
3305 let x: Vec<f64> = (0..60).map(|_| normal(&mut rng) * 0.4 - 0.3).collect();
3306 assert!(upper_bound(&x) < 0.0, "60 sessions of a clear harm trip it: {}", upper_bound(&x));
3307 }
3308
3309 fn monitored(free_helps: bool) -> Summary {
3313 summary::with_fit(monitored_sim(free_helps).0.summary(), None)
3314 }
3315
3316 fn monitored_sim(free_helps: bool) -> (Sim, Truth) {
3318 let mut t = truth();
3319 t.b[ASK] = 1.5;
3320 if free_helps {
3321 t.b[ASK] = 0.0;
3322 }
3323 let mut sim = Sim::new(31);
3324 let mut half = |grades: &[u8], _: bool, rng: &mut Rng| -> (Vec<Route>, Vec<f64>, Vec<f64>, &'static str) {
3325 let routes = grades.iter().map(|_| if rng.unit() < 0.5 { Route::Ask } else { Route::Free }).collect();
3326 (routes, vec![0.5; grades.len()], vec![0.5; grades.len()], "learned")
3327 };
3328 let mut sim_session = |sim: &mut Sim, name: &str| {
3330 for _ in 0..6 {
3331 let grades = vec![3, 2];
3332 let ids = vec![format!("L{}", sim.rng.next_u64() % 8), format!("L{}", sim.rng.next_u64() % 8)];
3333 sim.check(name, false, &grades, &ids, &t, &mut half);
3334 }
3335 };
3336 for i in 0..60 {
3337 sim_session(&mut sim, &format!("s{i}"));
3338 }
3339 (sim, t)
3340 }
3341
3342 #[test]
3348 fn a_placebo_neither_rule_delivers_moves_neither_the_estimate_nor_its_health() {
3349 const MEAN: f64 = -0.155_976_075_091_821;
3351 const UCB: f64 = -0.117_250_989_323_595;
3352 let (mut sim, t) = monitored_sim(false);
3353 let without = summary::with_fit(sim.summary(), None).monitor;
3354 let n = 12;
3357 let mut as_fixed = |grades: &[u8], _: bool, _: &mut Rng| -> (Vec<Route>, Vec<f64>, Vec<f64>, &'static str) {
3358 let routes: Vec<Route> = (0..grades.len())
3359 .map(|j| match j {
3360 j if j < conscience::CAP => Route::Free,
3361 j if j + 1 == grades.len() => Route::Placebo,
3362 _ => Route::Not,
3363 })
3364 .collect();
3365 let pf = routes.iter().map(|r| f64::from(u8::from(*r == Route::Free))).collect();
3366 (routes, vec![0.0; grades.len()], pf, "learned")
3367 };
3368 let ids: Vec<String> = (0..n).map(|i| format!("P{i}")).collect();
3369 sim.check("s59", false, &vec![1; n], &ids, &t, &mut as_fixed);
3370 let s = summary::with_fit(sim.summary(), None);
3371 let judged = |r: Route| s.applies.iter().filter(|a| a.route == r && matches!(a.label, summary::Label::Changed | summary::Label::Ignored)).count();
3372 assert_eq!(judged(Route::Placebo), 1, "the scenario: the placebo is judged");
3373 assert_eq!(fixed_routes(&ids.iter().map(|l| unit(l, 1, 0)).collect::<Vec<_>>(), &[]).iter().filter(|r| **r == Route::Free).count(), conscience::CAP);
3374 let m = &s.monitor;
3375 assert!((m.mean - MEAN).abs() < 1e-9, "mean {} vs 0.7.2's {MEAN}", m.mean);
3376 assert!(m.ucb.is_some_and(|u| (u - UCB).abs() < 1e-9), "ucb {:?} vs 0.7.2's {UCB}", m.ucb);
3377 assert_eq!((m.labelled, m.truncated), (without.labelled, without.truncated), "{m:?} vs {without:?}");
3379 assert!(m.max_weight < 10.0 && m.ess > 0.5 * f64::from(m.labelled), "{m:?}");
3380 }
3381
3382 #[test]
3387 fn one_label_does_not_carry_its_route() {
3388 for changed in [true, false] {
3389 let (mut sim, t) = monitored_sim(true);
3390 let without = summary::with_fit(sim.summary(), None).monitor;
3391 let mut deviated = |grades: &[u8], _: bool, _: &mut Rng| -> (Vec<Route>, Vec<f64>, Vec<f64>, &'static str) {
3392 (vec![Route::Ask, Route::Placebo], vec![1.0, 0.5], vec![0.0; grades.len()], "learned")
3393 };
3394 let truth = Truth { b: if changed { [9.0; P] } else { [-9.0; P] }, ..t };
3395 sim.check("s59", false, &[3, 2], &["Q0".to_string(), "Q1".to_string()], &truth, &mut deviated);
3396 let s = summary::with_fit(sim.summary(), None);
3397 let want = if changed { summary::Label::Changed } else { summary::Label::Ignored };
3398 assert_eq!(s.applies.iter().filter(|a| a.route == Route::Placebo).map(|a| a.label).collect::<Vec<_>>(), [want], "the scenario");
3399 let shift = (s.monitor.mean - without.mean).abs();
3400 assert!(shift < 0.05, "changed {changed}: one placebo moved the mean by {shift} ({} vs {})", s.monitor.mean, without.mean);
3401 }
3402 }
3403
3404 #[test]
3408 fn labels_without_leverage_do_not_open_a_route() {
3409 let n = 12;
3410 let mut as_fixed = |grades: &[u8], _: bool, _: &mut Rng| -> (Vec<Route>, Vec<f64>, Vec<f64>, &'static str) {
3411 let routes: Vec<Route> = (0..grades.len())
3412 .map(|j| match j {
3413 j if j < conscience::CAP => Route::Free,
3414 j if j + 1 == grades.len() => Route::Placebo,
3415 _ => Route::Not,
3416 })
3417 .collect();
3418 let pf = routes.iter().map(|r| f64::from(u8::from(*r == Route::Free))).collect();
3419 (routes, vec![0.0; grades.len()], pf, "learned")
3420 };
3421 let mut deviated = |grades: &[u8], _: bool, _: &mut Rng| -> (Vec<Route>, Vec<f64>, Vec<f64>, &'static str) {
3422 (vec![Route::Ask, Route::Placebo], vec![1.0, 0.5], vec![0.0; grades.len()], "learned")
3423 };
3424 for changed in [true, false] {
3425 let (mut sim, t) = monitored_sim(true);
3426 for k in 0..ROUTE_LABELS_MIN {
3427 let ids: Vec<String> = (0..n).map(|i| format!("Z{k}-{i}")).collect();
3428 sim.check(&format!("z{k}"), false, &vec![1; n], &ids, &t, &mut as_fixed);
3429 }
3430 let before = summary::with_fit(sim.summary(), None);
3431 let placebos = |s: &Summary| {
3432 s.applies.iter().filter(|a| a.route == Route::Placebo && matches!(a.label, summary::Label::Changed | summary::Label::Ignored)).count()
3433 };
3434 assert_eq!(placebos(&before), ROUTE_LABELS_MIN, "the scenario: every zero-leverage placebo judged");
3435 let truth = Truth { b: if changed { [9.0; P] } else { [-9.0; P] }, ..t };
3436 sim.check("s59", false, &[3, 2], &["Q0".to_string(), "Q1".to_string()], &truth, &mut deviated);
3437 let s = summary::with_fit(sim.summary(), None);
3438 assert_eq!(placebos(&s), ROUTE_LABELS_MIN + 1);
3439 let shift = (s.monitor.mean - before.monitor.mean).abs();
3440 assert!(shift < 0.05, "changed {changed}: one placebo moved the mean by {shift} ({} vs {})", s.monitor.mean, before.monitor.mean);
3441 }
3442 }
3443
3444 #[test]
3445 fn the_monitor_reverts_a_harmful_history_and_keeps_a_benign_one() {
3446 let harmful = monitored(false);
3447 assert!(harmful.monitor.reverted && harmful.monitor.ucb.is_some_and(|u| u < 0.0), "{:?}", harmful.monitor);
3448 assert_eq!(harmful.monitor.sessions, 60);
3449 let m = &harmful.monitor;
3452 assert!(m.labelled > 0 && m.labelled <= m.units, "{m:?}");
3453 assert!(m.ess > 0.0 && m.ess <= f64::from(m.labelled) + 1e-9 && m.max_weight >= 1.0 && m.truncated <= m.labelled, "{m:?}");
3454 let benign = monitored(true);
3455 assert!(!benign.monitor.reverted, "{:?}", benign.monitor);
3456 assert!(benign.monitor.mean > harmful.monitor.mean);
3457 }
3458
3459 #[test]
3462 fn the_fixed_routes_are_the_fixed_rule_by_position() {
3463 assert!(fixed_routes(&[], &[]).is_empty());
3464 let block = [unit("A", 1, 0), unit("B", 3, 1), unit("C", 2, 2)];
3465 assert_eq!(fixed_routes(&block, &[]), [Route::Ask; 3], "a 3 blocks: every 3, 2 and 1 up to the cap");
3466 let soft = [unit("A", 1, 0), unit("B", 2, 1)];
3467 assert_eq!(fixed_routes(&soft, &[]), [Route::Free; 2], "no 3: handed over free");
3468 let many: Vec<String> = (0..12).map(|i| format!("L{i}")).collect();
3469 let over: Vec<Unit> = many.iter().enumerate().map(|(r, l)| unit(l, 1, r)).collect();
3470 let routes = fixed_routes(&over, &[]);
3471 assert_eq!(routes.iter().filter(|r| **r == Route::Free).count(), conscience::CAP);
3472 assert_eq!(&routes[conscience::CAP..], &[Route::Not; 2], "past the cap: not now");
3473 let twice: Vec<Unit> = (0..11).map(|r| unit(if r < 9 { many[r].as_str() } else { "Z" }, 1, r)).collect();
3475 assert_eq!(&fixed_routes(&twice, &[])[9..], &[Route::Free, Route::Not]);
3476 }
3477
3478 #[test]
3479 fn a_prior_fit_writes_no_key() {
3480 assert_eq!(serde_json::to_string(&Fit::default()).unwrap(), "{}");
3481 assert_eq!(serde_json::from_str::<Fit>("{}").unwrap(), Fit::default());
3482 }
3483
3484 #[test]
3488 fn the_evidence_rule_is_chooses_without_the_test() {
3489 let f = free_beats_ask();
3490 let open = Gate { gain: vec![0.0; DRAWS], allowance: f64::INFINITY, ..Gate::default() };
3491 let ids: Vec<String> = (0..6).map(|i| format!("L{i}")).collect();
3492 let mut left = 0;
3493 for (n, g) in [(3, [2, 1, 1]), (3, [3, 2, 1]), (2, [1, 1, 0])] {
3494 let units: Vec<Unit> = (0..n).map(|r| unit(&ids[r], g[r], r)).collect();
3495 for stop in [false, true] {
3496 let check = Check { session: "new", stop, tokens: 0 };
3497 let routes = decide_on_evidence(&units, &check, &f, &[]);
3498 for seed in 0..5 {
3499 let d = choose(&units, &check, &f, &open, 0.0, seed, &[]);
3500 assert_eq!((&routes, d.xk), (&d.routes, None), "{g:?} stop {stop} seed {seed}");
3501 }
3502 left += usize::from(routes != fixed_routes(&units, &[]));
3503 assert_eq!(decide_on_evidence(&units, &check, &Fit::default(), &[]), fixed_routes(&units, &[]), "no model: the fixed rule");
3504 }
3505 }
3506 assert!(left > 0, "the scenario leaves the fixed rule");
3507 }
3508
3509 #[test]
3514 fn under_the_prior_the_evidence_rule_is_the_fixed_rule_on_every_grade_table() {
3515 let p = prior();
3516 let open = Gate { gain: vec![0.0; DRAWS], allowance: f64::INFINITY, ..Gate::default() };
3517 let ids: Vec<String> = (0..8).map(|i| format!("L{i}")).collect();
3518 let mut grades = tables(5);
3519 grades.extend([vec![2; 8], vec![1; 8], [vec![3], vec![2; 7]].concat()]);
3520 for g in &grades {
3521 let units: Vec<Unit> = ids.iter().zip(g).enumerate().map(|(r, (l, g))| unit(l, *g, r)).collect();
3522 for had in [Vec::new(), vec!["L1".to_string()]] {
3523 let want = fixed_routes(&units, &had);
3524 for stop in [true, false] {
3525 let check = Check { session: "s", stop, tokens: 0 };
3526 assert_eq!(decide_on_evidence(&units, &check, &p, &had), want, "{g:?} stop {stop} had {had:?}");
3527 let d = choose(&units, &check, &p, &open, 0.0, 7, &had);
3528 assert_eq!((&d.routes, d.source), (&want, if units.is_empty() { FIXED } else { EVIDENCE }), "{g:?} stop {stop} had {had:?}");
3529 }
3530 }
3531 }
3532 }
3533
3534 #[test]
3540 fn a_fit_with_evidence_of_dilution_deviates() {
3541 let world = |dilution: f64, sd: f64| {
3542 let mut beta = vec![0.0; P];
3543 beta[ASK] = 0.2;
3544 beta[NOT] = -1.0;
3545 beta[ASK_M] = dilution;
3546 let mut chol = vec![0.0; P * P];
3547 for k in 0..P {
3548 chol[k * P + k] = if k == ASK_M { sd } else { 0.01 };
3549 }
3550 let model = Model { beta, chol, sd: vec![1.0; P], spread: [0.01; FAMILIES], ..Model::default() };
3551 let rates = Rates { free_stop: 1.0, free_live: 1.0, later: 0.5, wait: 0.0, ask_live: 1.0 };
3552 Fit { key: "dilution".into(), model: Some(model), rates, ..Fit::default() }
3553 };
3554 let ids: Vec<String> = (0..9).map(|i| format!("L{i}")).collect();
3555 let units: Vec<Unit> = ids.iter().enumerate().map(|(r, l)| unit(l, if r == 0 { 3 } else { 2 }, r)).collect();
3556 let base = fixed_routes(&units, &[]);
3557 assert_eq!(base, [Route::Ask; 9], "the scenario: the fixed rule asks all nine");
3558 let gate = open_gate();
3559 for stop in [false, true] {
3560 let check = Check { session: "n", stop, tokens: 0 };
3561 let d = choose(&units, &check, &world(-1.0, 0.01), &gate, ALPHA, 1, &[]);
3562 assert_eq!((d.source, d.xk, d.xs), (EVIDENCE, None, None), "{d:?}");
3563 assert_ne!(d.routes, base, "{d:?}");
3564 assert!(d.routes.contains(&Route::Free) && d.routes.contains(&Route::Ask), "asks few, hands over the rest: {d:?}");
3565 assert_eq!(d.p, d.routes.iter().map(|r| f64::from(u8::from(*r == Route::Ask))).collect::<Vec<_>>(), "exact");
3566 assert_eq!(d.pf, d.routes.iter().map(|r| f64::from(u8::from(*r == Route::Free))).collect::<Vec<_>>(), "exact");
3567 let wide = world(-0.2, 0.5);
3569 let m = wide.model.as_ref().expect("a model");
3570 let feats = Feat::of(&units, stop, 0);
3571 let t = theta_mean(m, &feats, &units, &Levels::of(check.session, &m.judge, &units));
3572 let at_mode = solve(&Table::of(&Valuer { t: &t, feats: &feats, rates: &wide.rates }, &[false; 9]), crate::attention::RAIL, conscience::CAP);
3573 assert_ne!(at_mode, base, "stop {stop}");
3574 for f in [wide, world(0.0, 0.01)] {
3575 assert_eq!(decide_on_evidence(&units, &check, &f, &[]), base, "stop {stop}");
3576 }
3577 }
3578 }
3579
3580 fn tested_at(units: &[Unit], had: &[String], seeds: u64) -> Vec<Decision> {
3583 let p = prior();
3584 let open = Gate { gain: vec![0.0; DRAWS], allowance: 1e9, ..Gate::default() };
3585 (0..seeds).map(|seed| choose(units, &Check { session: "s", stop: false, tokens: 0 }, &p, &open, ALPHA, seed, had)).collect()
3586 }
3587
3588 #[test]
3593 fn the_tests_count_is_uniform_and_its_propensities_exact() {
3594 let ids: Vec<String> = (0..5).map(|i| format!("L{i}")).collect();
3595 let units: Vec<Unit> = ids.iter().zip([1, 2, 1, 2, 1]).enumerate().map(|(r, (l, g))| unit(l, g, r)).collect();
3596 let order = [1, 3, 0, 2, 4];
3597 let n = 1200;
3598 let mut counts = [0usize; 6];
3599 for d in tested_at(&units, &[], n) {
3600 let (Some(5), Some(s)) = (d.xk, d.xs) else { panic!("{d:?}") };
3601 let s = s as usize;
3602 counts[s] += 1;
3603 for (i, &j) in order.iter().enumerate() {
3604 assert_eq!(d.routes[j], if i < s { Route::Ask } else { Route::Free }, "{d:?}");
3605 assert_eq!((d.p[j], d.pf[j]), (outcome::round_p((5 - i) as f64 / 6.0), outcome::round_p((i + 1) as f64 / 6.0)), "{d:?}");
3606 }
3607 assert_eq!(d.source, if s == 0 { EXPERIMENT_CONTROL } else { EXPERIMENT }, "{d:?}");
3608 }
3609 let q = 1.0 / 6.0;
3610 let se = (q * (1.0 - q) / n as f64).sqrt();
3611 for (s, c) in counts.iter().enumerate() {
3612 let freq = *c as f64 / n as f64;
3613 assert!((freq - q).abs() <= 4.0 * se, "s = {s}: {freq:.4} against {q:.4} ({counts:?})");
3614 }
3615 }
3616
3617 #[test]
3620 fn the_test_never_touches_the_stop_or_a_three() {
3621 let ids: Vec<String> = (0..3).map(|i| format!("L{i}")).collect();
3622 let units: Vec<Unit> = ids.iter().zip([2, 3, 1]).enumerate().map(|(r, (l, g))| unit(l, g, r)).collect();
3623 let base = fixed_routes(&units, &[]);
3624 let p = prior();
3625 let open = Gate { gain: vec![0.0; DRAWS], allowance: 1e9, ..Gate::default() };
3626 for seed in 0..40 {
3627 let d = choose(&units, &Check { session: "s", stop: true, tokens: 0 }, &p, &open, ALPHA, seed, &[]);
3628 assert_eq!((&d.routes, d.source, d.xk, d.xs), (&base, EVIDENCE, None, None), "{d:?}");
3629 }
3630 let mut seen = std::collections::BTreeSet::new();
3631 for d in tested_at(&units, &[], 60) {
3632 assert_eq!((d.routes[1], d.p[1], d.pf[1], d.xk), (Route::Ask, 1.0, 0.0, Some(2)), "{d:?}");
3633 let s = d.xs.unwrap_or(9);
3634 assert_eq!(d.source == EXPERIMENT_CONTROL, s == 2, "{d:?}");
3635 seen.insert(s);
3636 }
3637 assert_eq!(seen.len(), 3, "every count drawn: {seen:?}");
3638 }
3639
3640 #[test]
3643 fn the_test_never_hands_over_a_lesson_the_session_had() {
3644 let ids: Vec<String> = (0..4).map(|i| format!("L{i}")).collect();
3645 let had = ["L2".to_string()];
3646 let soft: Vec<Unit> = ids.iter().zip([2, 1, 2, 1]).enumerate().map(|(r, (l, g))| unit(l, g, r)).collect();
3647 let block: Vec<Unit> = ids.iter().zip([2, 3, 2, 1]).enumerate().map(|(r, (l, g))| unit(l, g, r)).collect();
3648 for (units, kept) in [(&soft, Route::Not), (&block, Route::Ask)] {
3649 let mut moved = 0;
3650 for d in tested_at(units, &had, 60) {
3651 assert_eq!((d.routes[2], d.pf[2], d.xk), (kept, 0.0, Some(if kept == Route::Not { 3 } else { 2 })), "{d:?}");
3652 moved += usize::from(d.source == EXPERIMENT);
3653 }
3654 assert!(moved > 0, "the scenario: the test ran");
3655 }
3656 }
3657
3658 #[test]
3662 fn the_gate_refuses_the_test_when_short_and_risk_zero_turns_it_off() {
3663 let ids: Vec<String> = (0..3).map(|i| format!("L{i}")).collect();
3664 let units: Vec<Unit> = ids.iter().zip([2, 1, 1]).enumerate().map(|(r, (l, g))| unit(l, g, r)).collect();
3665 let base = fixed_routes(&units, &[]);
3666 let exact = |r: Route| base.iter().map(|x| f64::from(u8::from(*x == r))).collect::<Vec<_>>();
3667 let p = prior();
3668 let short = Gate { gain: vec![0.0; DRAWS], allowance: 1e-12, ..Gate::default() };
3669 let check = Check { session: "s", stop: false, tokens: 0 };
3670 let (mut refused, mut control) = (0, 0);
3671 for seed in 0..40 {
3672 let d = choose(&units, &check, &p, &short, ALPHA, seed, &[]);
3673 assert_eq!((&d.routes, &d.p, &d.pf, d.risk, d.xk), (&base, &exact(Route::Ask), &exact(Route::Free), 0.0, Some(3)), "{d:?}");
3674 match d.xs {
3675 Some(0) => control += usize::from(d.source == EXPERIMENT_CONTROL),
3676 _ => refused += usize::from(d.source == GATE),
3677 }
3678 let off = choose(&units, &check, &p, &open_gate(), 0.0, seed, &[]);
3679 assert_eq!((&off.routes, off.source, off.xk, off.xs), (&base, EVIDENCE, None, None), "{off:?}");
3680 }
3681 assert!(refused > 0 && control > 0 && refused + control == 40, "refused {refused}, control {control}");
3682 let room = Gate { allowance: 1e9 * ALPHA, ..short };
3685 assert!((0..40).any(|seed| choose(&units, &check, &p, &room, ALPHA, seed, &[]).source == EXPERIMENT));
3686 assert!((0..40).all(|seed| choose(&units, &check, &p, &room, 1e-12, seed, &[]).source != EXPERIMENT), "α scales it down");
3687 }
3688
3689 #[test]
3694 fn risk_zero_refuses_an_evidence_deviation_the_allowance_alone_would_cover() {
3695 let f = free_beats_ask();
3696 let ids: Vec<String> = (0..4).map(|i| format!("L{i}")).collect();
3698 let units: Vec<Unit> = ids.iter().enumerate().map(|(r, l)| unit(l, if r == 0 { 3 } else { 2 }, r)).collect();
3699 let base = fixed_routes(&units, &[]);
3700 let check = Check { session: "n", stop: false, tokens: 0 };
3701 let loss = 2.0 * units.len() as f64;
3704 let gate = Gate { gain: vec![-loss; DRAWS], allowance: 2.0 * loss, ..Gate::default() };
3705 let at = choose(&units, &check, &f, &gate, ALPHA, 1, &[]);
3706 assert_eq!(at.source, EVIDENCE, "the scenario: the allowance admits the deviation at ALPHA: {at:?}");
3707 assert_ne!(at.routes, base, "{at:?}");
3708 let off = choose(&units, &check, &f, &gate, 0.0, 1, &[]);
3709 assert_eq!((&off.routes, off.source, off.risk, off.gate, off.xk), (&base, GATE, 0.0, 0.0, None), "{off:?}");
3710 }
3711
3712 #[test]
3715 fn the_dilution_prior_has_its_mode_at_zero() {
3716 let c = coefs();
3717 for k in [ASK_M, FREE_M] {
3718 assert_eq!((c[k].mean, c[k].sd, c[k].nonpos, c[k].ard), (0.0, 0.1, true, false), "{}", NAMES[k]);
3719 }
3720 let p = prior();
3721 let m = p.model.as_ref().expect("the prior is a model");
3722 let units = [unit("A", 2, 0), unit("B", 2, 1), unit("C", 2, 2)];
3723 let feats = Feat::of(&units, false, 0);
3724 let t = theta_mean(m, &feats, &units, &Levels::of("s", &m.judge, &units));
3725 let v = Valuer { t: &t, feats: &feats, rates: &p.rates };
3726 assert!(v.ask(0, 1) > 0.0 && (v.ask(0, 1) - v.ask(0, 3)).abs() < 1e-12, "{} vs {}", v.ask(0, 1), v.ask(0, 3));
3727 assert!((v.free(0, 1) - v.free(0, 3)).abs() < 1e-12);
3728 }
3729
3730 #[test]
3734 fn a_bounded_slope_is_drawn_from_its_truncated_prior() {
3735 let p = prior();
3736 let m = p.model.as_ref().expect("the prior is a model");
3737 const N: u64 = 4000;
3738 let draws: Vec<[f64; P]> = (0..N).map(|k| draw_beta(m, gate_seed(&p, k as usize))).collect();
3739 let half = -DILUTION_PRIOR.1 * (2.0 / std::f64::consts::PI).sqrt();
3740 for k in [ASK_M, FREE_M] {
3742 assert!(draws.iter().all(|b| b[k] < 0.0), "{}: a draw at 0", NAMES[k]);
3743 let mean = draws.iter().map(|b| b[k]).sum::<f64>() / N as f64;
3744 assert!((mean - half).abs() < 0.005, "{}: mean {mean} against {half}", NAMES[k]);
3745 }
3746 let ask = draws.iter().map(|b| b[ASK]).sum::<f64>() / N as f64;
3747 assert!((ask - ASK_PRIOR.0).abs() < 0.03, "ask {ask}");
3748 for z in [-2.0, -0.3, 0.0, 1.0, 2.5] {
3750 assert!((below_zero(-1.0, 0.1, z) - (-1.0 + 0.1 * z)).abs() < 1e-6, "z {z}");
3751 let x = below_zero(50.0, 0.1, z);
3752 assert!(x <= 0.0 && x > -0.01, "z {z}: {x}");
3753 }
3754 assert!(below_zero(0.0, 0.1, -1.0) < below_zero(0.0, 0.1, 1.0), "monotone in z");
3755 assert!((phi(1.959_964) - 0.975).abs() < 1e-6 && (probit(0.975f64.ln(), 0.025) - 1.959_964).abs() < 1e-5);
3756 }
3757
3758 #[test]
3761 fn the_prior_is_a_model_of_the_priors_alone() {
3762 let p = prior();
3763 let m = p.model.as_ref().expect("the prior is a model");
3764 assert!(m.current() && m.converged && p.key.is_empty() && p.labels == 0);
3765 for (k, c) in coefs().iter().enumerate() {
3766 assert_eq!(m.beta[k], if c.nonpos { c.mean.min(0.0) } else { c.mean }, "{}", NAMES[k]);
3767 assert!((m.chol[k * P + k] - c.sd).abs() < 1e-12, "{}", NAMES[k]);
3768 }
3769 assert_eq!(m.spread, [LEVEL_START; FAMILIES]);
3770 assert_eq!(p.rates, Rates { free_stop: 0.5, free_live: 0.5, later: 0.5, wait: 0.5, ask_live: 1.0 });
3771 assert!(fit(&Summary::default(), None).model.is_none());
3772 }
3773
3774 #[test]
3777 fn a_prediction_is_the_posterior_chance_of_a_change() {
3778 let f = Feat { grade: 3, rank: 1, stop: true, three: true, pool: 0.0 };
3779 let mut beta = vec![0.0; P];
3780 beta[B0] = 0.8;
3781 beta[G3] = 0.4;
3782 beta[ASK] = 0.3;
3783 let lessons = BTreeMap::from([("A".to_string(), (0.5, 0.0))]);
3784 let certain = Model { beta: beta.clone(), chol: vec![0.0; P * P], sd: vec![1.0; P], spread: [0.0; FAMILIES], lessons, ..Model::default() };
3785 let mut fit = Fit { model: Some(certain.clone()), ..Fit::default() };
3786 let eta = dot(&std::array::from_fn(|k| beta[k]), &design(&f, Route::Ask, 2)) + 0.5;
3787 assert!((predict(&fit, &f, Route::Ask, 2, "s", "A").unwrap() - sigmoid(eta)).abs() < 1e-12);
3788 let mut wide = certain.clone();
3789 for k in 0..P {
3790 wide.chol[k * P + k] = 1.0;
3791 }
3792 wide.spread = [0.5; FAMILIES];
3793 fit.model = Some(wide);
3794 let p = predict(&fit, &f, Route::Ask, 2, "s", "unseen").unwrap();
3795 assert!(p > 0.5 && p < sigmoid(eta - 0.5), "{p} against {}", sigmoid(eta - 0.5));
3796 fit.model = Some(Model { beta: vec![0.0; P - 1], ..certain });
3797 assert_eq!(predict(&fit, &f, Route::Ask, 2, "s", "A"), None);
3798 }
3799}