pmcore/algorithms/
postprob.rs

1use crate::{
2    algorithms::Status,
3    prelude::algorithms::Algorithms,
4    structs::{
5        psi::{calculate_psi, Psi},
6        theta::Theta,
7    },
8};
9use anyhow::{Context, Result};
10use faer::Col;
11use pharmsol::prelude::{
12    data::{Data, ErrorModels},
13    simulator::Equation,
14};
15
16use crate::routines::evaluation::ipm::burke;
17use crate::routines::initialization;
18use crate::routines::output::{CycleLog, NPResult};
19use crate::routines::settings::Settings;
20
21/// Posterior probability algorithm
22/// Reweights the prior probabilities to the observed data and error model
23pub struct POSTPROB<E: Equation> {
24    equation: E,
25    psi: Psi,
26    theta: Theta,
27    w: Col<f64>,
28    objf: f64,
29    cycle: usize,
30    status: Status,
31    data: Data,
32    settings: Settings,
33    cyclelog: CycleLog,
34    error_models: ErrorModels,
35}
36
37impl<E: Equation> Algorithms<E> for POSTPROB<E> {
38    fn new(settings: Settings, equation: E, data: Data) -> Result<Box<Self>, anyhow::Error> {
39        Ok(Box::new(Self {
40            equation,
41            psi: Psi::new(),
42            theta: Theta::new(),
43            w: Col::zeros(0),
44            objf: f64::INFINITY,
45            cycle: 0,
46            status: Status::Starting,
47            error_models: settings.errormodels().clone(),
48            settings,
49            data,
50            cyclelog: CycleLog::new(),
51        }))
52    }
53    fn into_npresult(&self) -> NPResult<E> {
54        NPResult::new(
55            self.equation.clone(),
56            self.data.clone(),
57            self.theta.clone(),
58            self.psi.clone(),
59            self.w.clone(),
60            self.objf,
61            self.cycle,
62            self.status.clone(),
63            self.settings.clone(),
64            self.cyclelog.clone(),
65        )
66    }
67    fn get_settings(&self) -> &Settings {
68        &self.settings
69    }
70
71    fn equation(&self) -> &E {
72        &self.equation
73    }
74
75    fn get_data(&self) -> &Data {
76        &self.data
77    }
78
79    fn get_prior(&self) -> Theta {
80        initialization::sample_space(&self.settings).unwrap()
81    }
82
83    fn likelihood(&self) -> f64 {
84        self.objf
85    }
86
87    fn inc_cycle(&mut self) -> usize {
88        0
89    }
90
91    fn get_cycle(&self) -> usize {
92        0
93    }
94
95    fn set_theta(&mut self, theta: Theta) {
96        self.theta = theta;
97    }
98
99    fn theta(&self) -> &Theta {
100        &self.theta
101    }
102
103    fn psi(&self) -> &Psi {
104        &self.psi
105    }
106
107    fn set_status(&mut self, status: Status) {
108        self.status = status;
109    }
110
111    fn status(&self) -> &Status {
112        &self.status
113    }
114
115    fn convergence_evaluation(&mut self) {
116        // POSTPROB algorithm converges after a single evaluation
117        self.status = Status::MaxCycles;
118    }
119
120    fn converged(&self) -> bool {
121        true
122    }
123
124    fn evaluation(&mut self) -> Result<()> {
125        self.psi = calculate_psi(
126            &self.equation,
127            &self.data,
128            &self.theta,
129            &self.error_models,
130            false,
131            false,
132        )?;
133        (self.w, self.objf) = burke(&self.psi).context("Error in IPM")?;
134        Ok(())
135    }
136
137    fn condensation(&mut self) -> Result<()> {
138        Ok(())
139    }
140    fn optimizations(&mut self) -> Result<()> {
141        Ok(())
142    }
143
144    fn logs(&self) {}
145
146    fn expansion(&mut self) -> Result<()> {
147        Ok(())
148    }
149}