pmcore/estimation/nonparametric/
summaries.rs

1use ndarray::{Array1, Array2};
2use pharmsol::{Data, Equation, Event};
3
4use crate::estimation::nonparametric::NonParametricResult;
5use crate::estimation::nonparametric::{population_mean_median, posterior_mean_median};
6use crate::results::{FitSummary, IndividualSummary, ParameterSummary, PopulationSummary};
7
8pub fn fit_summary<E: Equation>(result: &NonParametricResult<E>) -> FitSummary {
9    FitSummary {
10        objective_function: result.objf(),
11        converged: result.converged(),
12        iterations: result.cycles(),
13        subject_count: result.data().subjects().len(),
14        observation_count: count_observations(result.data()),
15        parameter_count: result.get_theta().parameters().len(),
16    }
17}
18
19pub fn population_summary<E: Equation>(result: &NonParametricResult<E>) -> PopulationSummary {
20    let theta_matrix = to_ndarray_matrix(result.get_theta().matrix());
21    let weights = Array1::from_iter(result.weights().iter());
22    let (mean, median) = population_mean_median(&theta_matrix, &weights)
23        .expect("population summary should be derivable from theta and weights");
24
25    let parameters = result
26        .get_theta()
27        .parameters()
28        .names()
29        .into_iter()
30        .enumerate()
31        .map(|(index, name)| {
32            let column = theta_matrix.column(index).to_vec();
33            let mean_value = mean[index];
34            let sd = weighted_sd(&column, &weights, mean_value);
35            let cv_percent = if mean_value.abs() > f64::EPSILON {
36                (sd / mean_value.abs()) * 100.0
37            } else {
38                0.0
39            };
40
41            ParameterSummary {
42                name,
43                mean: mean_value,
44                median: median[index],
45                sd,
46                cv_percent,
47            }
48        })
49        .collect();
50
51    PopulationSummary { parameters }
52}
53
54pub fn individual_summaries<E: Equation>(
55    result: &NonParametricResult<E>,
56) -> Vec<IndividualSummary> {
57    let theta_matrix = to_ndarray_matrix(result.get_theta().matrix());
58    let psi_matrix = to_ndarray_matrix(result.psi().matrix());
59    let weights = Array1::from_iter(result.weights().iter());
60    let (means, _) = posterior_mean_median(&theta_matrix, &psi_matrix, &weights)
61        .expect("individual summaries should be derivable from theta, psi, and weights");
62    let parameter_names = result.get_theta().parameters().names();
63
64    result
65        .data()
66        .subjects()
67        .iter()
68        .enumerate()
69        .map(|(subject_index, subject)| IndividualSummary {
70            id: subject.id().clone(),
71            parameter_names: parameter_names.clone(),
72            estimates: means.row(subject_index).to_vec(),
73            standard_errors: None,
74        })
75        .collect()
76}
77
78fn count_observations(data: &Data) -> usize {
79    data.subjects()
80        .iter()
81        .flat_map(|subject| subject.occasions())
82        .flat_map(|occasion| occasion.events())
83        .filter(|event| matches!(event, Event::Observation(_)))
84        .count()
85}
86
87fn to_ndarray_matrix(matrix: &faer::Mat<f64>) -> Array2<f64> {
88    Array2::from_shape_fn((matrix.nrows(), matrix.ncols()), |(row, col)| {
89        matrix[(row, col)]
90    })
91}
92
93fn weighted_sd(values: &[f64], weights: &Array1<f64>, mean: f64) -> f64 {
94    let variance = values
95        .iter()
96        .zip(weights.iter())
97        .map(|(value, weight)| weight * (value - mean).powi(2))
98        .sum::<f64>();
99    variance.sqrt()
100}
pmcore/estimation/nonparametric/summaries.rs

pmcore/estimation/nonparametric/
summaries.rs
