34 lines
1.2 KiB
Rust
34 lines
1.2 KiB
Rust
/// calculate the binomial coefficient.
|
|
pub fn nchoosek(n: u64, k: u64) -> f64 {
|
|
if k > n {
|
|
return 0.0;
|
|
}
|
|
let n_fac = rug::Integer::from(rug::Integer::factorial(n as u32)) / (rug::Integer::from(rug::Integer::factorial(k as u32)) * rug::Integer::from(rug::Integer::factorial((n - k) as u32)));
|
|
return n_fac.to_f64();
|
|
}
|
|
|
|
/// the probability that a given outcome (with probability p) will occur at least k times in n independent trials (aka Bernoulli trials).
|
|
pub fn at_least_with_replacement(k: u64, n: u64, p: f64) -> rug::Float {
|
|
let mut prob_at_least = rug::Float::with_val(64, 0.0);
|
|
for x in k..n {
|
|
prob_at_least += nchoosek(n, x) * rug::Float::with_val(64, p.powi(x as i32)) * rug::Float::with_val(64, (1.0 - p).powi((n - x) as i32));
|
|
}
|
|
prob_at_least
|
|
}
|
|
|
|
#[cfg(test)]
|
|
mod tests {
|
|
use crate::probability::binomial::{at_least_with_replacement, nchoosek};
|
|
|
|
#[test]
|
|
fn test_nchoosek() {
|
|
assert_eq!(17310309456440u64, nchoosek(100, 10));
|
|
}
|
|
|
|
#[test]
|
|
fn test_at_least() {
|
|
// probability of at least 48 heads in 100 coin flips
|
|
assert_eq!(0.691350293205374, at_least_with_replacement(48, 100, 0.5).to_f64());
|
|
}
|
|
}
|