diff --git a/src/estimation/nonparametric/theta.rs b/src/estimation/nonparametric/theta.rs index f0504c677..615a4bdd4 100644 --- a/src/estimation/nonparametric/theta.rs +++ b/src/estimation/nonparametric/theta.rs @@ -112,11 +112,16 @@ impl Theta { } /// Suggest a new support point to add to the matrix - /// The point is only added if it is at least `min_dist` away from all existing support points - /// and within the limits specified by `limits` + /// The point is first clamped to the parameter bounds, and is only added if it is + /// at least `min_dist` away from all existing support points pub(crate) fn suggest_point(&mut self, spp: &[f64], min_dist: f64) -> Result<()> { - if self.check_point(spp, min_dist) { - self.add_point(spp)?; + let spp: Vec = spp + .iter() + .zip(self.parameters.finite_ranges()) + .map(|(val, (lower, upper))| val.clamp(lower, upper)) + .collect(); + if self.check_point(&spp, min_dist) { + self.add_point(&spp)?; } Ok(()) } @@ -525,6 +530,19 @@ mod tests { .add("v", 5.0, 50.0) } + #[test] + fn suggest_point_clamps_to_bounds() { + let params = parameters(); + let mut theta = + Theta::from_parts(Mat::from_fn(1, 2, |_, j| [0.5, 20.0][j]), params).unwrap(); + + theta.suggest_point(&[6.8, 1.0], 1e-4).unwrap(); + + assert_eq!(theta.nspp(), 2); + assert_eq!(*theta.matrix().get(1, 0), 1.0); + assert_eq!(*theta.matrix().get(1, 1), 5.0); + } + fn temp_csv_path() -> String { format!("test_temp_theta_{}.csv", rand::random::()) } diff --git a/tests/onecomp.rs b/tests/onecomp.rs index c63509b52..d74fdc00d 100644 --- a/tests/onecomp.rs +++ b/tests/onecomp.rs @@ -131,9 +131,17 @@ fn test_one_compartment_npod() -> Result<()> { .fit_with(NonParametricAlgorithm::npod())?; // Check the results - assert_eq!(result.cycles(), 11); + assert_eq!(result.cycles(), 6); assert!(result.objf() - 565.7749 < 0.01); + // Support points stay inside the parameter box + let theta = result.get_theta().matrix(); + for (j, (lower, upper)) in parameters.finite_ranges().into_iter().enumerate() { + for i in 0..theta.nrows() { + assert!((lower..=upper).contains(&theta[(i, j)])); + } + } + Ok(()) }