Support points as defined structures

benches/performance.rs

@@ -38,8 +38,9 @@ fn readme(n: usize) {
         },
         (2, 1),
     );
+    let spp = support_point!["ka" => 0.3, "ke" => 0.5, "tlag" => 0.1, "v" => 70.0];
     for _ in 0..n {
-        let op = ode.estimate_predictions(&subject, &vec![0.3, 0.5, 0.1, 70.0]);
+        let op = ode.estimate_predictions(&subject, &spp);
         assert_eq!(
             op.flat_predictions(),
             vec![

examples/bke.rs

@@ -19,7 +19,7 @@ fn main() {
         |_p| fa! {},
         |_p, _t, _cov, _x| {},
         |x, p, _t, _cov, y| {
-            fetch_params!(p, _ke, v);
+            let v = p.get("v").unwrap();
             y[0] = x[0] / v;
         },
         (1, 1),
@@ -28,39 +28,31 @@ fn main() {
     let ode = equation::ODE::new(
         |x, p, _t, dx, rateiv, _cov| {
             // fetch_cov!(cov, t, wt);
-            fetch_params!(p, ke, _v);
+            let ke = p.get("ke").unwrap();
             dx[0] = -ke * x[0] + rateiv[0];
         },
         |_p| lag! {},
         |_p| fa! {},
         |_p, _t, _cov, _x| {},
         |x, p, _t, _cov, y| {
-            fetch_params!(p, _ke, v);
+            let v = p.get("v").unwrap();
             y[0] = x[0] / v;
         },
         (1, 1),
     );
 
-    let em = ErrorModel::new((0.0, 0.05, 0.0, 0.0), 0.0, &ErrorType::Add);
-    let ll = an.estimate_likelihood(
-        &subject,
-        &vec![1.02282724609375, 194.51904296875],
-        &em,
-        false,
-    );
-    let op = an.estimate_predictions(&subject, &vec![1.02282724609375, 194.51904296875]);
+    let em = ErrorModel::new((0.0, 0.05, 0.0, 0.0), 0.0, &ErrorType::Additive);
+    // let spp = SupportPoint::new(vec![("ke", 1.02282724609375), ("v", 194.51904296875)]);
+    let spp = support_point!("ke" => 1.02282724609375, "v" => 194.51904296875);
+    let ll = an.estimate_likelihood(&subject, &spp, &em, false);
+    let op = an.estimate_predictions(&subject, &spp);
     println!(
         "Analytical: \n-2ll:{:#?}\n{:#?}",
         -2.0 * ll,
         op.flat_predictions()
     );
 
-    let ll = ode.estimate_likelihood(
-        &subject,
-        &vec![1.02282724609375, 194.51904296875],
-        &em,
-        false,
-    );
-    let op = ode.estimate_predictions(&subject, &vec![1.02282724609375, 194.51904296875]);
+    let ll = ode.estimate_likelihood(&subject, &spp, &em, false);
+    let op = ode.estimate_predictions(&subject, &spp);
     println!("ODE: \n-2ll:{:#?}\n{:#?}", -2.0 * ll, op.flat_predictions());
 }

examples/error.rs

@@ -26,7 +26,7 @@ fn main() {
         |_p, _t, _cov, _x| {},
         |x, p, t, cov, y| {
             fetch_cov!(cov, t, WT);
-            fetch_params!(p, CL0, V0, Vp0, Q0);
+            let V0 = p.get("v0").unwrap();
             let V = V0 / (WT / 85.0);
             y[0] = x[0] / V;
         },
@@ -39,6 +39,9 @@ fn main() {
         .repeat(120, 0.1)
         .build();
 
-    let op = ode.estimate_predictions(&subject, &vec![0.1, 0.1, 0.1, 0.1, 70.0]);
+    //spp = vec![0.1, 0.1, 0.1, 0.1, 70.0]
+    let spp = support_point!("cl0" => 0.1, "v0" => 0.1, "vp0" => 0.1, "q0" => 0.1, "v0" => 70.0);
+
+    let op = ode.estimate_predictions(&subject, &spp);
     dbg!(op);
 }

examples/gendata.rs

@@ -23,8 +23,7 @@ fn main() {
             // user defined
             dx[0] = -ke * x[0];
         },
-        |p, d| {
-            fetch_params!(p, _ke0);
+        |_p, d| {
             d[1] = 0.1;
         },
         |_p| lag! {},
@@ -33,8 +32,7 @@ fn main() {
             fetch_params!(p, ke0);
             x[1] = ke0;
         },
-        |x, p, _t, _cov, y| {
-            fetch_params!(p, _ke0);
+        |x, _p, _t, _cov, y| {
             y[0] = x[0] / 50.0;
         },
         (2, 1),
@@ -62,6 +60,7 @@ fn main() {
 
     let mut data = vec![];
     for (i, spp) in support_points.iter().enumerate() {
+        let spp = SupportPoint::from_vec(spp.clone(), vec!["ke", "ske", "v"]);
         let trajectories = sde.estimate_predictions(&subject, &spp);
         let trajectory = trajectories.row(0);
         // dbg!(&trajectory);

examples/nsim.rs

@@ -1,5 +1,6 @@
 fn main() {
     use pharmsol::*;
+
     let subject = Subject::builder("id1")
         .bolus(0.0, 100.0, 0)
         .repeat(2, 0.5)
@@ -15,7 +16,7 @@ fn main() {
     let ode = equation::ODE::new(
         |x, p, t, dx, _rateiv, cov| {
             fetch_cov!(cov, t, wt, age);
-            fetch_params!(p, ka, ke, _tlag, _v);
+            fetch_params!(p, ka, ke);
             // Secondary Eqs
 
             let ke = ke * wt.powf(0.75) * (age / 25.0).powf(0.5);
@@ -25,13 +26,13 @@ fn main() {
             dx[1] = ka * x[0] - ke * x[1];
         },
         |p| {
-            fetch_params!(p, _ka, _ke, tlag, _v);
+            fetch_params!(p, tlag);
             lag! {0=>tlag}
         },
         |_p| fa! {},
         |_p, _t, _cov, _x| {},
         |x, p, _t, _cov, y| {
-            fetch_params!(p, _ka, _ke, _tlag, v);
+            fetch_params!(p, v);
             y[0] = x[1] / v;
         },
         (2, 1),
@@ -61,6 +62,13 @@ fn main() {
     //     (2, 1),
     // );
 
-    let op = ode.estimate_predictions(&subject, &vec![0.3, 0.5, 0.1, 70.0]);
+    let spp = support_point!(
+        "ka" => 0.3,
+        "ke" => 0.5,
+        "tlag" => 0.1,
+        "v" => 70.0
+    );
+
+    let op = ode.estimate_predictions(&subject, &spp);
     println!("{op:#?}");
 }

examples/ode_readme.rs

@@ -17,24 +17,26 @@ fn main() {
     let ode = equation::ODE::new(
         |x, p, _t, dx, _rateiv, _cov| {
             // fetch_cov!(cov, t,);
-            fetch_params!(p, ka, ke, _tlag, _v);
+            fetch_params!(p, ka, ke);
             //Struct
             dx[0] = -ka * x[0];
             dx[1] = ka * x[0] - ke * x[1];
         },
         |p| {
-            fetch_params!(p, _ka, _ke, tlag, _v);
+            fetch_params!(p, tlag);
             lag! {0=>tlag}
         },
         |_p| fa! {},
         |_p, _t, _cov, _x| {},
         |x, p, _t, _cov, y| {
-            fetch_params!(p, _ka, _ke, _tlag, v);
+            fetch_params!(p, v);
             y[0] = x[1] / v;
         },
         (2, 1),
     );
 
-    let op = ode.estimate_predictions(&subject, &vec![0.3, 0.5, 0.1, 70.0]);
+    let spp = support_point!["ka" => 0.3, "ke" => 0.5, "tlag" => 0.1, "v" => 70.0];
+
+    let op = ode.estimate_predictions(&subject, &spp);
     println!("{:#?}", op.flat_predictions());
 }

examples/pf.rs

@@ -12,8 +12,9 @@ fn main() {
 
     let sde = equation::SDE::new(
         |x, p, _t, dx, _rateiv, _cov| {
+            fetch_params!(p, ke0);
             dx[0] = -x[0] * x[1]; // ke *x[0]
-            dx[1] = -x[1] + p[0]; // mean reverting
+            dx[1] = -x[1] + ke0; // mean reverting
         },
         |_p, d| {
             d[0] = 1.0;
@@ -29,10 +30,12 @@ fn main() {
         10000,
     );
 
-    let em = ErrorModel::new((0.5, 0.0, 0.0, 0.0), 0.0, &error_model::ErrorType::Add); // sigma = 0.5
+    let em = ErrorModel::new((0.5, 0.0, 0.0, 0.0), 0.0, &error_model::ErrorType::Additive); // sigma = 0.5
 
-    let ll = sde.estimate_likelihood(&subject, &vec![1.0], &em, false);
+    let spp = support_point!["ke0" => 1.0];
 
-    dbg!(sde.estimate_likelihood(&subject, &vec![1.0], &em, false));
+    let ll = sde.estimate_likelihood(&subject, &spp, &em, false);
+
+    dbg!(sde.estimate_likelihood(&subject, &spp, &em, false));
     println!("{ll:#?}");
 }

examples/sim_theta.rs

@@ -0,0 +1,61 @@
+use csv::ReaderBuilder;
+use pharmsol::*;
+use prelude::{data::read_pmetrics, simulator::Prediction};
+
+fn main() {
+    // path to theta
+    let path = "../PMcore/outputs/theta.csv";
+    // read theta into an Array2<f64>
+    let mut rdr = ReaderBuilder::new().from_path(path).unwrap();
+
+    let mut spps = vec![];
+    for result in rdr.records() {
+        let record = result.unwrap();
+        let mut row = vec![];
+        for field in record.iter() {
+            row.push(field.parse::<f64>().unwrap());
+        }
+        spps.push(row);
+    }
+
+    let sde = equation::SDE::new(
+        |x, p, _t, dx, _rateiv, _cov| {
+            // automatically defined
+            fetch_params!(p, ke0);
+            // let ke0 = 1.2;
+            dx[1] = -x[1] + ke0;
+            let ke = x[1];
+            // user defined
+            dx[0] = -ke * x[0];
+        },
+        |p, d| {
+            fetch_params!(p, _ke0);
+            d[1] = 0.1;
+        },
+        |_p| lag! {},
+        |_p| fa! {},
+        |p, _t, _cov, x| {
+            fetch_params!(p, ke0);
+            x[1] = ke0;
+        },
+        |x, p, _t, _cov, y| {
+            fetch_params!(p, _ke0);
+            y[0] = x[0] / 50.0;
+        },
+        (2, 1),
+        1,
+    );
+
+    let data = read_pmetrics("../PMcore/examples/iov/test.csv").unwrap();
+
+    for (i, spp) in spps.iter().enumerate() {
+        for (j, subject) in data.get_subjects().iter().enumerate() {
+            let spp = SupportPoint::from_vec(spp.clone(), vec!["ke0"]);
+            let trajectories: ndarray::Array2<Prediction> =
+                sde.estimate_predictions(&subject, &spp);
+            let trajectory = trajectories.row(0);
+            println!("{}, {}", i, j);
+            dbg!(trajectory);
+        }
+    }
+}

src/data/covariate.rs

@@ -1,14 +1,14 @@
-use serde::Deserialize;
+use serde::{Deserialize, Serialize};
 use std::{collections::HashMap, fmt};
 
-#[derive(Clone, Debug, Deserialize)]
+#[derive(Clone, Debug, Deserialize, Serialize)]
 pub enum InterpolationMethod {
     Linear { slope: f64, intercept: f64 },
     CarryForward { value: f64 },
 }
 
 /// A [CovariateSegment] is a segment of the piece-wise interpolation of a [Covariate]
-#[derive(Clone, Debug, Deserialize)]
+#[derive(Clone, Debug, Deserialize, Serialize)]
 pub struct CovariateSegment {
     from: f64,
     to: f64,
@@ -46,7 +46,7 @@ impl CovariateSegment {
 }
 
 /// A [Covariate] is a collection of [CovariateSegment]s, which allows for interpolation of covariate values
-#[derive(Clone, Debug, Deserialize)]
+#[derive(Clone, Debug, Deserialize, Serialize)]
 pub struct Covariate {
     name: String,
     segments: Vec<CovariateSegment>,
@@ -113,7 +113,7 @@ impl fmt::Display for Covariate {
 }
 
 /// [Covariates] is a collection of [Covariate]
-#[derive(Clone, Debug, Deserialize)]
+#[derive(Clone, Debug, Deserialize, Serialize)]
 pub struct Covariates {
     // Mapping from covariate name to its segments
     // FIXME: this hashmap have a key, covariate also has a name field, we are never

src/data/error_model.rs

@@ -1,3 +1,5 @@
+use serde::{Deserialize, Serialize};
+
 use crate::simulator::likelihood::Prediction;
 
 #[derive(Debug, Clone)]
@@ -16,10 +18,10 @@ impl<'a> ErrorModel<'a> {
     }
 }
 
-#[derive(Debug, Clone)]
+#[derive(Debug, Clone, Serialize, Deserialize, Copy)]
 pub enum ErrorType {
-    Add,
-    Prop,
+    Additive,
+    Proportional,
 }
 #[allow(clippy::extra_unused_lifetimes)]
 impl<'a> ErrorModel<'_> {
@@ -34,8 +36,8 @@ impl<'a> ErrorModel<'_> {
             + c3 * prediction.observation().powi(3);
 
         let res = match self.e_type {
-            ErrorType::Add => (alpha.powi(2) + self.gl.powi(2)).sqrt(),
-            ErrorType::Prop => self.gl * alpha,
+            ErrorType::Additive => (alpha.powi(2) + self.gl.powi(2)).sqrt(),
+            ErrorType::Proportional => self.gl * alpha,
         };
 
         if res.is_nan() || res < 0.0 {

src/data/event.rs

@@ -1,9 +1,9 @@
 use std::fmt;
 
-use serde::Deserialize;
+use serde::{Deserialize, Serialize};
 
 /// An Event can be a Bolus, Infusion, or Observation
-#[derive(Debug, Clone, Deserialize)]
+#[derive(Debug, Clone, Deserialize, Serialize)]
 pub enum Event {
     Bolus(Bolus),
     Infusion(Infusion),
@@ -29,7 +29,7 @@ impl Event {
 /// An instantaenous input of drug
 ///
 /// An amount of drug is added to a compartment at a specific time
-#[derive(Debug, Clone, Deserialize)]
+#[derive(Debug, Clone, Deserialize, Serialize)]
 pub struct Bolus {
     time: f64,
     amount: f64,
@@ -62,7 +62,7 @@ impl Bolus {
 }
 
 /// A continuous dose of drug
-#[derive(Debug, Clone, Deserialize)]
+#[derive(Debug, Clone, Deserialize, Serialize)]
 pub struct Infusion {
     time: f64,
     amount: f64,
@@ -99,7 +99,7 @@ impl Infusion {
 }
 
 /// An observation of drug concentration or covariates
-#[derive(Debug, Clone, Deserialize)]
+#[derive(Debug, Clone, Deserialize, Serialize)]
 pub struct Observation {
     time: f64,
     value: f64,