Merge pull request #5 from aseyboldt/performance

Some more performance improvements

src/lib.rs

@@ -200,22 +200,18 @@ pub fn run_with(
         panic!("x and y must have at least one element");
     }
 
-    let x0 = x[0];
+    // Use the middle of the time period as reference
+    // to prevent strong correlations between alpha and beta
+    let x0 = x.iter().sum::<f64>() / x.len() as f64;
 
     let x = x.iter().map(|x| x - x0).collect::<Vec<_>>();
 
     let model = Regression::new(x.clone(), y.clone());
 
     // y = alpha + beta * x + noise
-    let guessed_beta = y.iter().sum::<f64>() / x.iter().sum::<f64>();
+    let guessed_beta = 0.;//y.iter().sum::<f64>() / x.iter().sum::<f64>();
     let guessed_alpha = y.iter().sum::<f64>() / y.len() as f64;
-    let guessed_sigma = x
-        .iter()
-        .zip(y.iter())
-        .map(|(x, y)| (y - guessed_alpha - guessed_beta * x).powi(2))
-        .sum::<f64>()
-        .sqrt()
-        / y.len() as f64;
+    let guessed_sigma = 1.;
     let initial_position = vec![guessed_alpha, guessed_beta, guessed_sigma];
     log(format!("initial_position = {:?}", initial_position).as_str());
 

src/model/regression.rs

@@ -91,8 +91,8 @@ impl CpuLogpFunc for Regression {
         let logp_beta = log_pdf_normal_propto(beta, 10f64.ln(), 0.01);
         let logp_sigma = 0.; // flat prior
 
-        let mut d_logp_d_alpha = 0.;
-        let mut d_logp_d_beta = 0.;
+        let mut d_logp_d_alpha = -alpha / 100.;
+        let mut d_logp_d_beta = -beta / 100.;
         let mut d_logp_d_sigma = 0.;
 
         let mut logp_y = 0.;

src/sampler.rs

@@ -59,7 +59,7 @@ where
         .expect("Unrecoverable error during init");
 
     // Burn the first x samples to get away from the initial position
-    for _ in 0..50 {
+    for _ in 0..num_tune {
         sampler.draw().expect("Unrecoverable error during burning");
     }