Skip to main content

pmcore/estimation/nonparametric/sampling/
latin.rs

1use anyhow::Result;
2use faer::Mat;
3use rand::prelude::*;
4use rand::rngs::StdRng;
5
6use crate::estimation::nonparametric::Theta;
7use crate::model::{BoundedParameter, ParameterSpace};
8
9pub fn generate(
10    parameters: &ParameterSpace<BoundedParameter>,
11    points: usize,
12    seed: usize,
13) -> Result<Theta> {
14    let ranges = parameters.finite_ranges();
15    let mut rng = StdRng::seed_from_u64(seed as u64);
16
17    let mut intervals = Vec::new();
18    for _ in 0..ranges.len() {
19        let mut param_intervals: Vec<f64> = (0..points).map(|i| i as f64).collect();
20        param_intervals.shuffle(&mut rng);
21        intervals.push(param_intervals);
22    }
23
24    let rand_matrix = Mat::from_fn(points, ranges.len(), |i, j| {
25        let interval = intervals[j][i];
26        let random_offset = rng.random::<f64>();
27        let unscaled = (interval + random_offset) / points as f64;
28        let (lower, upper) = ranges[j];
29        lower + unscaled * (upper - lower)
30    });
31
32    Theta::from_parts(rand_matrix, parameters.clone())
33}
34
35#[cfg(test)]
36mod tests {
37    use super::*;
38    use crate::model::{BoundedParameter, ParameterSpace};
39
40    #[test]
41    fn latin_generate_produces_requested_shape() {
42        let params = ParameterSpace::<BoundedParameter>::new()
43            .add("a", 0.0, 1.0)
44            .add("b", 0.0, 1.0)
45            .add("c", 0.0, 1.0);
46
47        let theta = generate(&params, 10, 22).unwrap();
48        assert_eq!(theta.nspp(), 10);
49        assert_eq!(theta.matrix().ncols(), 3);
50    }
51}