Skip to content
Merged
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
67 changes: 50 additions & 17 deletions src/distribution/internal.rs
Original file line number Diff line number Diff line change
Expand Up @@ -258,13 +258,18 @@
use crate::distribution::{Continuous, ContinuousCDF, Discrete, DiscreteCDF};
use crate::prec;

enum CheckContinuousError {
LnPdf { x: f64 },
PdfStep { x: f64, _sum: f64 },
}

/// cdf should be the integral of the pdf
fn check_integrate_pdf_is_cdf<D: ContinuousCDF<f64, f64> + Continuous<f64, f64>>(
dist: &D,
x_min: f64,
x_max: f64,
step: f64,
) {
) -> Result<f64, CheckContinuousError> {
let mut prev_x = x_min;
let mut prev_density = dist.pdf(x_min);
let mut sum = 0.0;
Expand All @@ -276,19 +281,18 @@
assert!(density >= 0.0);

let ln_density = dist.ln_pdf(x);
if ln_density.is_finite()
&& !prec::abs_diff_eq!(density.ln(), ln_density, epsilon = 1e-10)
{
return Err(CheckContinuousError::LnPdf { x });

Check warning on line 287 in src/distribution/internal.rs

View check run for this annotation

Codecov / codecov/patch

src/distribution/internal.rs#L287

Added line #L287 was not covered by tests
}

prec::assert_abs_diff_eq!(density.ln(), ln_density, epsilon = 1e-10);

// triangle rule
// trapezoidal rule
sum += (prev_density + density) * step / 2.0;

let cdf = dist.cdf(x);
if !prec::abs_diff_eq!(sum, cdf, epsilon = 1e-3) {
panic!(
"Integral of pdf doesn't equal cdf!\n\
Integration from {x_min} by {step} to {x} = {sum}\n\
cdf = {cdf}"
);
return Err(CheckContinuousError::PdfStep { x, _sum: sum });
}

if x >= x_max {
Expand All @@ -299,10 +303,7 @@
}
}

assert!(
sum > 0.99 && sum <= 1.001,
"sum should be close to 1, but is {sum}"
);
Ok(sum)
}

/// cdf should be the sum of the pmf
Expand Down Expand Up @@ -341,7 +342,7 @@
x_min: f64,
x_max: f64,
step: f64,
) {
) -> Result<(), f64> {
const DELTA: f64 = 1e-12;
const DX: f64 = 2.0 * DELTA;
let mut prev_x = x_min;
Expand All @@ -354,14 +355,18 @@

let d_cdf = dist.cdf(x_ahead) - dist.cdf(x_behind);

prec::assert_abs_diff_eq!(d_cdf, DX * density, epsilon = 1e-11);
if !prec::abs_diff_eq!(d_cdf, DX * density, epsilon = 1e-11) {
return Err(x);

Check warning on line 359 in src/distribution/internal.rs

View check run for this annotation

Codecov / codecov/patch

src/distribution/internal.rs#L359

Added line #L359 was not covered by tests
}

if x >= x_max {
break;
} else {
prev_x = x;
}
}

Ok(())
}

/// Does a series of checks that all continuous distributions must obey.
Expand All @@ -379,8 +384,36 @@
assert_eq!(dist.cdf(f64::NEG_INFINITY), 0.0);
assert_eq!(dist.cdf(f64::INFINITY), 1.0);

check_integrate_pdf_is_cdf(dist, x_min, x_max, (x_max - x_min) / 100000.0);
check_derivative_of_cdf_is_pdf(dist, x_min, x_max, (x_max - x_min) / 100000.0);
let integrate_res =
check_integrate_pdf_is_cdf(dist, x_min, x_max, (x_max - x_min) / 100000.0);
let diff_res =
check_derivative_of_cdf_is_pdf(dist, x_min, x_max, (x_max - x_min) / 100000.0);
match integrate_res {
// if integration failed along the way...
Err(CheckContinuousError::LnPdf { x })

Check warning on line 393 in src/distribution/internal.rs

View check run for this annotation

Codecov / codecov/patch

src/distribution/internal.rs#L393

Added line #L393 was not covered by tests
| Err(CheckContinuousError::PdfStep { x, .. }) => {
// then check if differentiation fails along the way
if let Err(diff_err_x) = diff_res {
panic!(
"integration mismatched around x = {x} \n\
and derivative insufficiently close at x={diff_err_x}"

Check warning on line 399 in src/distribution/internal.rs

View check run for this annotation

Codecov / codecov/patch

src/distribution/internal.rs#L397-L399

Added lines #L397 - L399 were not covered by tests
)
}
// otherwise: passes differentiation test
}
// if integration is close but not at endpoint (surprising!)...
Ok(sum) if !prec::abs_diff_eq!(sum, 1.0, epsilon = 1e-3) => {
// then check if differentiation fails along the way
if let Err(diff_err_x) = diff_res {
panic!(
"integration summed to {sum}, insufficiently close to 1.0 \n\
and derivative insufficiently close at x={diff_err_x}"

Check warning on line 410 in src/distribution/internal.rs

View check run for this annotation

Codecov / codecov/patch

src/distribution/internal.rs#L408-L410

Added lines #L408 - L410 were not covered by tests
)
}
}
// passes integration test
_ => (),
}
}

/// Does a series of checks that all positive discrete distributions must
Expand Down
Loading