Skip to content
Merged
3 changes: 0 additions & 3 deletions benches/vonmises.rs
Original file line number Diff line number Diff line change
@@ -1,17 +1,14 @@
use criterion::black_box;
use criterion::measurement::WallTime;
use criterion::AxisScale;
use criterion::BatchSize;
use criterion::BenchmarkId;
use criterion::Criterion;
use criterion::PlotConfiguration;
use criterion::Throughput;
use criterion::{criterion_group, criterion_main};
use rand::Rng;
use rv::dist::VonMises;
use rv::misc::bessel::log_i0;
use rv::prelude::*;
use rv::traits::*;
use std::f64::consts::PI;

fn bench_vm_draw(c: &mut Criterion) {
Expand Down
19 changes: 7 additions & 12 deletions examples/betaprime_sbc.rs
Original file line number Diff line number Diff line change
@@ -1,19 +1,14 @@
use rv::dist::BetaPrime;

#[cfg(feature = "experimental")]
use rand::SeedableRng;
#[cfg(feature = "experimental")]
use rand_xoshiro::Xoshiro256Plus;
#[cfg(feature = "experimental")]
use rv::experimental::stick_breaking_process::{
StickBreaking, StickBreakingDiscrete, StickBreakingDiscreteSuffStat,
};
use rv::prelude::*;

// Simulation-based calibration
// For details see http://www.stat.columbia.edu/~gelman/research/unpublished/sbc.pdf
#[cfg(feature = "experimental")]
fn main() {
use rand::SeedableRng;
use rand_xoshiro::Xoshiro256Plus;
use rv::experimental::stick_breaking_process::{
StickBreaking, StickBreakingDiscrete, StickBreakingDiscreteSuffStat,
};
use rv::prelude::*;

let mut rng = Xoshiro256Plus::seed_from_u64(123);
let n_samples = 10000;
let n_obs = 10;
Expand Down
14 changes: 6 additions & 8 deletions examples/sbd.rs
Original file line number Diff line number Diff line change
@@ -1,14 +1,12 @@
use rand::SeedableRng;
use rv::prelude::*;

#[cfg(feature = "experimental")]
use rv::experimental::stick_breaking_process::{
StickBreaking, StickBreakingDiscrete, StickSequence,
};

fn main() {
#[cfg(feature = "experimental")]
{
use rand_xoshiro::rand_core::SeedableRng;
use rv::experimental::stick_breaking_process::{
StickBreaking, StickBreakingDiscrete, StickSequence,
};
use rv::prelude::*;

// Instantiate a stick-breaking process
let alpha = 10.0;
let sbp = StickBreaking::new(UnitPowerLaw::new(alpha).unwrap());
Expand Down
12 changes: 5 additions & 7 deletions examples/stickbreaking_posterior.rs
Original file line number Diff line number Diff line change
@@ -1,13 +1,11 @@
use itertools::Either;
use peroxide::statistics::stat::Statistics;
use rv::prelude::*;

#[cfg(feature = "experimental")]
use rv::experimental::stick_breaking_process::*;

fn main() {
#[cfg(feature = "experimental")]
{
use itertools::Either;
use peroxide::fuga::Statistics;
use rv::experimental::stick_breaking_process::*;
use rv::prelude::*;

let mut rng = rand::thread_rng();
let sb = StickBreaking::new(UnitPowerLaw::new(3.0).unwrap());

Expand Down
7 changes: 3 additions & 4 deletions src/data/mod.rs
Original file line number Diff line number Diff line change
Expand Up @@ -235,7 +235,7 @@ pub fn extract_stat_then<X, Fx, Pr, Fnx, Y>(
) -> Y
where
Fx: HasSuffStat<X> + HasDensity<X>,
Pr: ConjugatePrior<X, Fx>,
Pr: ConjugatePrior<X, Fx> + ?Sized,
Fnx: Fn(&Fx::Stat) -> Y,
{
match x {
Expand Down Expand Up @@ -386,8 +386,8 @@ mod tests {
fn impl_bool_into_bool() {
let t = true;
let f = false;
assert_eq!(t.into_bool(), true);
assert_eq!(f.into_bool(), false);
assert!(t.into_bool());
assert!(!f.into_bool());
}

#[test]
Expand Down Expand Up @@ -578,7 +578,6 @@ mod tests {
&self,
_x: &crate::data::GaussianData<f64>,
) -> Self::PpCache {
()
}

fn ln_pp_with_cache(
Expand Down
8 changes: 4 additions & 4 deletions src/dist/betaprime.rs
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,9 @@ use std::f64;
use std::fmt;
use std::sync::OnceLock;

#[cfg(feature = "experimental")]
use super::UnitPowerLaw;

/// [Beta prime distribution](https://en.wikipedia.org/wiki/Beta_prime_distribution),
/// BetaPrime(α, β) over x in (0, ∞).
///
Expand Down Expand Up @@ -277,18 +280,14 @@ impl Sampleable<f64> for BetaPrime {
}
}

use crate::data::DataOrSuffStat;
#[cfg(feature = "experimental")]
use crate::experimental::stick_breaking_process::{
StickBreakingDiscrete, StickBreakingDiscreteSuffStat,
};
use crate::traits::ConjugatePrior;

#[cfg(feature = "experimental")]
use crate::experimental::stick_breaking_process::StickBreaking;

use crate::prelude::UnitPowerLaw;

#[cfg(feature = "experimental")]
impl Sampleable<StickBreakingDiscrete> for BetaPrime {
fn draw<R: Rng>(&self, rng: &mut R) -> StickBreakingDiscrete {
Expand Down Expand Up @@ -390,6 +389,7 @@ impl ConjugatePrior<usize, StickBreakingDiscrete> for BetaPrime {
) -> Self {
match data {
DataOrSuffStat::Data(xs) => {
#[allow(clippy::useless_asref)]
let stat = StickBreakingDiscreteSuffStat::from(xs.as_ref());
self.posterior(&DataOrSuffStat::SuffStat(&stat))
}
Expand Down
11 changes: 7 additions & 4 deletions src/dist/cdvm.rs
Original file line number Diff line number Diff line change
Expand Up @@ -243,10 +243,13 @@ impl HasSuffStat<usize> for Cdvm {
// TODO: Should we cache twopimu_over_m.cos() and twopimu_over_m.sin()?

let (sin_twopimu_over_m, cos_twopimu_over_m) = twopimu_over_m.sin_cos();
self.kappa
* (stat.sum_cos() * cos_twopimu_over_m
+ stat.sum_sin() * sin_twopimu_over_m)
- stat.n() as f64 * self.log_norm_const()
self.kappa.mul_add(
stat.sum_cos().mul_add(
cos_twopimu_over_m,
stat.sum_sin() * sin_twopimu_over_m,
),
-(stat.n() as f64 * self.log_norm_const()),
)
}
}

Expand Down
4 changes: 4 additions & 0 deletions src/dist/mod.rs
Original file line number Diff line number Diff line change
Expand Up @@ -45,7 +45,9 @@ mod pareto;
mod poisson;
mod scaled;
mod scaled_inv_chi_squared;
mod scaled_prior;
mod shifted;
mod shifted_prior;
mod skellam;
mod students_t;
mod uniform;
Expand Down Expand Up @@ -108,7 +110,9 @@ pub use scaled_inv_chi_squared::{
ScaledInvChiSquared, ScaledInvChiSquaredError,
ScaledInvChiSquaredParameters,
};
pub use scaled_prior::{ScaledPrior, ScaledPriorError};
pub use shifted::{Shifted, ShiftedError, ShiftedParameters};
pub use shifted_prior::{ShiftedPrior, ShiftedPriorError};
pub use skellam::{Skellam, SkellamError, SkellamParameters};
pub use students_t::{StudentsT, StudentsTError, StudentsTParameters};
pub use uniform::{Uniform, UniformError, UniformParameters};
Expand Down
30 changes: 15 additions & 15 deletions src/dist/normal_gamma/gaussian_prior.rs
Original file line number Diff line number Diff line change
Expand Up @@ -3,7 +3,7 @@ use std::f64::consts::LN_2;

use super::dos_to_post;
use crate::consts::*;
use crate::data::{extract_stat, GaussianSuffStat};
use crate::data::{extract_stat_then, GaussianSuffStat};
use crate::dist::{Gaussian, NormalGamma};
use crate::gaussian_prior_geweke_testable;
use crate::misc::ln_gammafn;
Expand Down Expand Up @@ -77,20 +77,20 @@ impl ConjugatePrior<f64, Gaussian> for NormalGamma {
}

fn ln_pp_cache(&self, x: &DataOrSuffStat<f64, Gaussian>) -> Self::PpCache {
let stat = extract_stat(self, x);

let params = posterior_from_stat(self, &stat);
let PosteriorParameters { r, s, v, .. } = params;

let half_v = v / 2.0;
let g_ratio = ln_gammafn(half_v + 0.5) - ln_gammafn(half_v);
let term = 0.5_f64.mul_add(LN_2, -HALF_LN_2PI)
+ 0.5_f64.mul_add(
(r / (r + 1_f64)).ln(),
half_v.mul_add(s.ln(), g_ratio),
);

(params, term)
extract_stat_then(self, x, |stat| {
let params = posterior_from_stat(self, stat);
let PosteriorParameters { r, s, v, .. } = params;

let half_v = v / 2.0;
let g_ratio = ln_gammafn(half_v + 0.5) - ln_gammafn(half_v);
let term = 0.5_f64.mul_add(LN_2, -HALF_LN_2PI)
+ 0.5_f64.mul_add(
(r / (r + 1_f64)).ln(),
half_v.mul_add(s.ln(), g_ratio),
);

(params, term)
})
}

fn ln_pp_with_cache(&self, cache: &Self::PpCache, y: &f64) -> f64 {
Expand Down
34 changes: 18 additions & 16 deletions src/dist/normal_inv_chi_squared/gaussian_prior.rs
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@ use std::collections::BTreeMap;
use std::f64::consts::PI;

use crate::consts::HALF_LN_PI;
use crate::data::{extract_stat, extract_stat_then, GaussianSuffStat};
use crate::data::{extract_stat_then, GaussianSuffStat};
use crate::dist::{Gaussian, NormalInvChiSquared};
use crate::gaussian_prior_geweke_testable;
use crate::misc::ln_gammafn;
Expand Down Expand Up @@ -69,10 +69,11 @@ impl ConjugatePrior<f64, Gaussian> for NormalInvChiSquared {
GaussianSuffStat::new()
}

fn posterior(&self, x: &DataOrSuffStat<f64, Gaussian>) -> Self {
extract_stat_then(self, x, |stat: &GaussianSuffStat| {
posterior_from_stat(self, &stat).into()
})
fn posterior_from_suffstat(
&self,
stat: &GaussianSuffStat,
) -> Self::Posterior {
posterior_from_stat(self, stat).into()
}

#[inline]
Expand All @@ -87,23 +88,24 @@ impl ConjugatePrior<f64, Gaussian> for NormalInvChiSquared {
) -> f64 {
extract_stat_then(self, x, |stat: &GaussianSuffStat| {
let n = stat.n() as f64;
let post: Self = posterior_from_stat(self, &stat).into();
let post: Self = posterior_from_stat(self, stat).into();
let lnz_n = post.ln_z();
n.mul_add(-HALF_LN_PI, lnz_n - cache)
})
}

fn ln_pp_cache(&self, x: &DataOrSuffStat<f64, Gaussian>) -> Self::PpCache {
let stat = extract_stat(self, x);
let post = posterior_from_stat(self, &stat);
let kn = post.kn;
let vn = post.vn;

let z = 0.5_f64.mul_add(
(kn / ((kn + 1.0) * PI * vn * post.s2n)).ln(),
ln_gammafn((vn + 1.0) / 2.0) - ln_gammafn(vn / 2.0),
);
(post, z)
extract_stat_then(self, x, |stat: &GaussianSuffStat| {
let post = posterior_from_stat(self, stat);
let kn = post.kn;
let vn = post.vn;

let z = 0.5_f64.mul_add(
(kn / ((kn + 1.0) * PI * vn * post.s2n)).ln(),
ln_gammafn((vn + 1.0) / 2.0) - ln_gammafn(vn / 2.0),
);
(post, z)
})
}

fn ln_pp_with_cache(&self, cache: &Self::PpCache, y: &f64) -> f64 {
Expand Down
Loading