From 86f822bf34fcd2caaf24d40196ba16f3429f61da Mon Sep 17 00:00:00 2001 From: Felipe Contreras Salinas Date: Mon, 8 Dec 2025 15:49:41 -0300 Subject: [PATCH] binomial p2 --- locales/en/main.ftl | 2 +- locales/es/main.ftl | 2 +- src/calc.rs | 418 ++++++++++++++++++++++++++--------- src/components/calculator.rs | 9 +- 4 files changed, 316 insertions(+), 115 deletions(-) diff --git a/locales/en/main.ftl b/locales/en/main.ftl index 5a77f1b..684bfbc 100644 --- a/locales/en/main.ftl +++ b/locales/en/main.ftl @@ -10,4 +10,4 @@ hyper-description = Hypergeometric distribution measures the probability of gett success-probability = Success probability trials-number = Number of trials successes-number = Number of successes -binom-description = . +binom-description = Binomial distribution measures the probability of getting a given amount of successes in a sequence of experiments. For example, if you flip 5 (number of trials = 5) balanced coins (success probability = 0.5), the distribution describes the probability of having a given number of heads (successes number = X). diff --git a/locales/es/main.ftl b/locales/es/main.ftl index a7ceb83..67928cb 100644 --- a/locales/es/main.ftl +++ b/locales/es/main.ftl @@ -10,4 +10,4 @@ hyper-description = La distribución hipergeométrica mide la probabilidad de ob success-probability = Probabilidad de éxito trials-number = Número de intentos successes-number = Número de éxitos -binom-description = . +binom-description = La distribución binomial mide la probabilidad de obtener un cierto número de éxitos en una secuencia de experimentos. Por ejemplo, si lanzas 5 (número de intentos = 5) monedas justas (probabilidad de éxito = 0.5), la distribución describe la probabilidad de obtener un cierto número de monedas (número de éxitos = X). diff --git a/src/calc.rs b/src/calc.rs index af21bb0..29829ed 100644 --- a/src/calc.rs +++ b/src/calc.rs @@ -2,79 +2,179 @@ use std::{collections::HashMap, iter::repeat}; -#[derive(Default)] +#[derive(Debug)] +pub struct HyperGeometricInput { + population_size: u8, + successes: u8, + sample_size: u8, + sample_successes: u8, +} + +impl HyperGeometricInput { + pub fn new( + population_size: u8, + successes: u8, + sample_size: u8, + sample_successes: u8, + ) -> Option { + if successes > population_size + || sample_size > population_size + || sample_successes > sample_size + { + None + } else { + Some(Self { + population_size, + successes, + sample_size, + sample_successes, + }) + } + } +} + +/// Result of hypergeometric probability calculation. +#[derive(Default, Debug, PartialEq)] pub struct HyperGeometricProb { + /// Probability of getting exactly X successes in the sample. pub exactly: f64, + /// Probability of getting strictly less than X successes in the sample. pub less_than: f64, + /// Probability of getting less than or exactly X successes in the sample. pub less_or_equal: f64, + /// Probability of getting strictly more X successes in the sample. pub greater_than: f64, + /// Probability of getting more than or exactly X successes in the sample. pub greater_or_equal: f64, } -pub fn hyper_geometric( - population_size: u8, - successes: u8, - sample_size: u8, - sample_successes: u8, -) -> Option { - if successes > population_size - || sample_size > population_size - || sample_successes > sample_size +pub fn hyper_geometric(input: HyperGeometricInput) -> HyperGeometricProb { + let exactly = hyper_geometric_exactly(&input); + let (less_than, less_or_equal, greater_or_equal, greater_than) = if input.sample_successes + < input.sample_size / 2 { - None + let less_than = (0..input.sample_successes) + .map(|i| { + hyper_geometric_exactly(&HyperGeometricInput { + population_size: input.population_size, + successes: input.successes, + sample_size: input.sample_size, + sample_successes: i, + }) + }) + .sum::() + .abs(); + let less_or_equal = less_than + exactly; + let greater_or_equal = (1.0 - less_than).abs(); + let greater_than = (1.0 - less_or_equal).abs(); + (less_than, less_or_equal, greater_or_equal, greater_than) } else { - let exactly = - hyper_geometric_exactly(population_size, successes, sample_size, sample_successes); - let (less_than, less_or_equal, greater_or_equal, greater_than) = - if sample_successes < sample_size / 2 { - let less_than = (0..sample_successes) - .map(|i| hyper_geometric_exactly(population_size, successes, sample_size, i)) - .sum::() - .abs(); - let less_or_equal = less_than + exactly; - let greater_or_equal = (1.0 - less_than).abs(); - let greater_than = (1.0 - less_or_equal).abs(); - (less_than, less_or_equal, greater_or_equal, greater_than) - } else { - let greater_than = (sample_successes + 1..=sample_size) - .map(|i| hyper_geometric_exactly(population_size, successes, sample_size, i)) - .sum::() - .abs(); - let greater_or_equal = greater_than + exactly; - let less_or_equal = (1.0 - greater_than).abs(); - let less_than = (1.0 - greater_or_equal).abs(); - (less_than, less_or_equal, greater_or_equal, greater_than) - }; - Some(HyperGeometricProb { - exactly, - less_than, - less_or_equal, - greater_than, - greater_or_equal, - }) + let greater_than = (input.sample_successes + 1..=input.sample_size.min(input.successes)) + .map(|i| { + hyper_geometric_exactly(&HyperGeometricInput { + population_size: input.population_size, + successes: input.successes, + sample_size: input.sample_size, + sample_successes: i, + }) + }) + .sum::() + .abs(); + let greater_or_equal = greater_than + exactly; + let less_or_equal = (1.0 - greater_than).abs(); + let less_than = (1.0 - greater_or_equal).abs(); + (less_than, less_or_equal, greater_or_equal, greater_than) + }; + HyperGeometricProb { + exactly, + less_than, + less_or_equal, + greater_than, + greater_or_equal, } } -#[derive(Default)] -pub struct BinomialProb { - pub exactly: f64, - pub less_than: f64, - pub less_or_equal: f64, - pub greater_than: f64, - pub greater_or_equal: f64, -} - -pub fn binomial( +pub struct BinomialInput { success_probability: f64, trials_number: u8, successes_number: u8, -) -> Option { - if successes_number > trials_number { - None - } else { - let (p_powers, pc_powers) = powers(success_probability, trials_number); - let exactly = binom_exactly(success_probability, trials_number, successes_number); - None +} + +impl BinomialInput { + pub fn new(success_probability: f64, trials_number: u8, successes_number: u8) -> Option { + if success_probability < 0.0 + || success_probability > 1.0 + || successes_number > trials_number + { + None + } else { + Some(Self { + success_probability, + trials_number, + successes_number, + }) + } + } +} + +#[derive(Default, Debug, PartialEq)] +pub struct BinomialProb { + pub exactly: f64, + pub less_than: f64, + pub less_or_equal: f64, + pub greater_than: f64, + pub greater_or_equal: f64, +} + +pub fn binomial(input: BinomialInput) -> BinomialProb { + let (p_powers, pc_powers) = powers(input.success_probability, input.trials_number); + let exactly = binomial_exactly(&input, &p_powers, &pc_powers); + let (less_than, less_or_equal, greater_or_equal, greater_than) = + if input.successes_number < input.trials_number / 2 { + let less_than = (0..input.successes_number) + .map(|i| { + binomial_exactly( + &BinomialInput { + success_probability: input.success_probability, + trials_number: input.trials_number, + successes_number: i, + }, + &p_powers, + &pc_powers, + ) + }) + .sum::() + .abs(); + let less_or_equal = less_than + exactly; + let greater_or_equal = (1.0 - less_than).abs(); + let greater_than = (1.0 - less_or_equal).abs(); + (less_than, less_or_equal, greater_or_equal, greater_than) + } else { + let greater_than = (input.successes_number + 1..=input.trials_number) + .map(|i| { + binomial_exactly( + &BinomialInput { + success_probability: input.success_probability, + trials_number: input.trials_number, + successes_number: i, + }, + &p_powers, + &pc_powers, + ) + }) + .sum::() + .abs(); + let greater_or_equal = greater_than + exactly; + let less_or_equal = (1.0 - greater_than).abs(); + let less_than = (1.0 - greater_or_equal).abs(); + (less_than, less_or_equal, greater_or_equal, greater_than) + }; + BinomialProb { + exactly, + less_than, + less_or_equal, + greater_than, + greater_or_equal, } } @@ -83,65 +183,49 @@ pub fn binomial( /// /// The formula is choose(successes, sample_successes) * choose(population_size - successes, /// sample_size - sample_successes) / choose(population_size, sample_size) -fn hyper_geometric_exactly( - population_size: u8, - successes: u8, - sample_size: u8, - sample_successes: u8, -) -> f64 { - if population_size == successes { - return if sample_successes == sample_size { +fn hyper_geometric_exactly(input: &HyperGeometricInput) -> f64 { + if input.population_size == input.successes { + return if input.sample_successes == input.sample_size { 1.0 } else { 0.0 }; } - if successes == 0 { - return if sample_successes == 0 { 1.0 } else { 0.0 }; + if input.successes == 0 { + return if input.sample_successes == 0 { + 1.0 + } else { + 0.0 + }; } // On top we have: successes!, (population_size - successes)!, sample_size! and // (population_size - sample_size)! - let top_factors = (1..=successes) - .chain(1..=(population_size - successes)) - .chain(1..=sample_size) - .chain(1..=(population_size - sample_size)) + let top_factors = (1..=input.successes) + .chain(1..=(input.population_size - input.successes)) + .chain(1..=input.sample_size) + .chain(1..=(input.population_size - input.sample_size)) .flat_map(|n| factorize(n)) - .fold(HashMap::::new(), |mut counts, i| { - *counts.entry(i).or_default() += 1; - counts - }); + .fold(HashMap::::new(), group_factors); // On bottom we have: sample_successes!, (successes - sample_successes)! // (sample_size - sample_successes)!, (population_size - successes - sample_size + sample_successes)! // and population_size! - let bot_factors = (1..=sample_successes) - .chain(1..=(successes - sample_successes)) - .chain(1..=(sample_size - sample_successes)) + let bot_factors = (1..=input.sample_successes) + .chain(1..=(input.successes - input.sample_successes)) + .chain(1..=(input.sample_size - input.sample_successes)) .chain( - 1..=((population_size as u16 + sample_successes as u16 - - successes as u16 - - sample_size as u16) as u8), + 1..=((input.population_size as u16 + input.sample_successes as u16 + - input.successes as u16 + - input.sample_size as u16) as u8), ) - .chain(1..=population_size) + .chain(1..=input.population_size) .flat_map(|n| factorize(n)) - .fold(HashMap::::new(), |mut counts, i| { - counts.entry(i).and_modify(|count| *count += 1).or_insert(1); - counts - }); + .fold(HashMap::::new(), group_factors); let (top_factors, bot_factors) = simplify(top_factors, bot_factors); - let top_product: f64 = top_factors - .into_iter() - .flat_map(|(f, count)| repeat(f).take(count as usize)) - .map(|f| f as f64) - .product(); - - let bot_product: f64 = bot_factors - .into_iter() - .flat_map(|(f, count)| repeat(f).take(count as usize)) - .map(|f| f as f64) - .product(); + let top_product = product(top_factors); + let bot_product = product(bot_factors); top_product / bot_product } @@ -149,19 +233,55 @@ fn hyper_geometric_exactly( /// Computes the probability of getting exactly `successes_number` within `trials_number` given /// that the success probability is `success_probability`. /// -/// The formula is choose(successes_number, trials_number) * (success_probability)^successes_number +/// The formula is choose(trials_number, successes_number) * (success_probability)^successes_number /// * (1 - success_probability)^(trials_number - successes_number) -fn binom_exactly(success_probability: f64, trials_number: u8, successes_number: u8) -> f64 { - 0.0 +fn binomial_exactly(input: &BinomialInput, p_powers: &[f64], pc_powers: &[f64]) -> f64 { + if input.success_probability == 0.0 { + return if input.successes_number == 0 { + 1.0 + } else { + 0.0 + }; + } + if input.success_probability == 1.0 { + return if input.successes_number == input.trials_number { + 1.0 + } else { + 0.0 + }; + } + choose(input.trials_number, input.successes_number) + * p_powers[input.successes_number as usize] + * pc_powers[(input.trials_number - input.successes_number) as usize] } -fn powers(p: f64, N: u8) -> (Vec, Vec) { - let mut p_powers = Vec::with_capacity((N + 1) as usize); - let mut pc_powers = Vec::with_capacity((N + 1) as usize); +fn choose(n: u8, k: u8) -> f64 { + // On top we have: n! + let top_factors = (1..=n) + .flat_map(|n| factorize(n)) + .fold(HashMap::::new(), group_factors); + + // On bottom we have: k!, (n - k)! + let bot_factors = (1..=k) + .chain(1..=(n - k)) + .flat_map(|n| factorize(n)) + .fold(HashMap::::new(), group_factors); + + let (top_factors, bot_factors) = simplify(top_factors, bot_factors); + + let top_product = product(top_factors); + let bot_product = product(bot_factors); + + top_product / bot_product +} + +fn powers(p: f64, n: u8) -> (Vec, Vec) { + let mut p_powers = Vec::with_capacity((n + 1) as usize); + let mut pc_powers = Vec::with_capacity((n + 1) as usize); let mut p_power = 1.0; let mut pc_power = 1.0; - for _ in 0..N + 1 { + for _ in 0..n + 1 { p_powers.push(p_power); pc_powers.push(pc_power); p_power = p_power * p; @@ -170,6 +290,19 @@ fn powers(p: f64, N: u8) -> (Vec, Vec) { (p_powers, pc_powers) } +fn group_factors(mut counts: HashMap, i: u8) -> HashMap { + *counts.entry(i).or_default() += 1; + counts +} + +fn product(factors: HashMap) -> f64 { + factors + .into_iter() + .flat_map(|(f, count)| repeat(f).take(count as usize)) + .map(|f| f as f64) + .product() +} + /// Simplify factors for a fraction. /// /// This assumes factors are already prime factors. @@ -259,7 +392,10 @@ fn factorize(n: u8) -> FactorIter<'static> { #[cfg(test)] mod test { - use crate::calc::hyper_geometric_exactly; + use crate::calc::{ + BinomialInput, BinomialProb, HyperGeometricInput, HyperGeometricProb, binomial, + binomial_exactly, hyper_geometric, hyper_geometric_exactly, powers, + }; use super::factorize; @@ -280,18 +416,80 @@ mod test { #[test] fn test_hypergeometric_exact_all_successes() { - assert_eq!(hyper_geometric_exactly(10, 10, 5, 5), 1.0); - assert_eq!(hyper_geometric_exactly(10, 10, 5, 4), 0.0); + let input = &HyperGeometricInput::new(10, 10, 5, 5).unwrap(); + assert_eq!(hyper_geometric_exactly(input), 1.0); + let input = &HyperGeometricInput::new(10, 10, 5, 4).unwrap(); + assert_eq!(hyper_geometric_exactly(input), 0.0); } #[test] fn test_hypergeometric_exact_no_successes() { - assert_eq!(hyper_geometric_exactly(10, 0, 5, 0), 1.0); - assert_eq!(hyper_geometric_exactly(10, 0, 5, 1), 0.0); + let input = &HyperGeometricInput::new(10, 0, 5, 0).unwrap(); + assert_eq!(hyper_geometric_exactly(input), 1.0); + let input = &HyperGeometricInput::new(10, 0, 5, 1).unwrap(); + assert_eq!(hyper_geometric_exactly(input), 0.0); } #[test] fn test_hypergeometric_exact() { - assert_eq!(hyper_geometric_exactly(10, 3, 5, 2), 5.0 / 12.0); + let input = &HyperGeometricInput::new(10, 3, 5, 2).unwrap(); + assert_eq!(hyper_geometric_exactly(input), 5.0 / 12.0); + } + + #[test] + fn test_hypergeometric_aces_poker() { + let input = HyperGeometricInput::new(52, 4, 5, 4).unwrap(); + let exact = 1.846892603195124e-5; + assert_eq!( + hyper_geometric(input), + HyperGeometricProb { + exactly: exact, + less_than: 1.0 - exact, + less_or_equal: 1.0, + greater_than: 0.0, + greater_or_equal: exact + } + ); + } + + #[test] + fn test_binom_exact_all_success() { + let (p_powers, pc_powers) = powers(1.0, 5); + let input = &BinomialInput::new(1.0, 5, 5).unwrap(); + assert_eq!(binomial_exactly(input, &p_powers, &pc_powers), 1.0); + let input = &BinomialInput::new(1.0, 5, 4).unwrap(); + assert_eq!(binomial_exactly(input, &p_powers, &pc_powers), 0.0); + } + + #[test] + fn test_binom_exact_no_success() { + let (p_powers, pc_powers) = powers(0.0, 5); + let input = &BinomialInput::new(0.0, 5, 0).unwrap(); + assert_eq!(binomial_exactly(input, &p_powers, &pc_powers), 1.0); + let input = &BinomialInput::new(0.0, 5, 1).unwrap(); + assert_eq!(binomial_exactly(input, &p_powers, &pc_powers), 0.0); + } + + #[test] + fn test_binomial_exact() { + let (p_powers, pc_powers) = powers(0.5, 5); + let input = &BinomialInput::new(0.5, 5, 3).unwrap(); + assert_eq!(binomial_exactly(input, &p_powers, &pc_powers), 10.0 / 32.0); + } + + #[test] + fn test_binomial() { + // 10.0 / 32.0 + let input = BinomialInput::new(0.5, 5, 3).unwrap(); + assert_eq!( + binomial(input), + BinomialProb { + exactly: 10.0 / 32.0, + less_than: 16.0 / 32.0, + less_or_equal: 26.0 / 32.0, + greater_than: 6.0 / 32.0, + greater_or_equal: 16.0 / 32.0, + } + ); } } diff --git a/src/components/calculator.rs b/src/components/calculator.rs index ac104a5..db9793c 100644 --- a/src/components/calculator.rs +++ b/src/components/calculator.rs @@ -6,7 +6,7 @@ use leptos::prelude::{ use leptos::{IntoView, component, view}; use leptos_fluent::move_tr; -use crate::calc::{binomial, hyper_geometric}; +use crate::calc::{BinomialInput, HyperGeometricInput, binomial, hyper_geometric}; #[component] pub fn HyperCalculator() -> impl IntoView { @@ -15,12 +15,13 @@ pub fn HyperCalculator() -> impl IntoView { let (sample, set_sample) = signal(0u8); let (sample_successes, set_sample_successes) = signal(0u8); let result = move || { - hyper_geometric( + HyperGeometricInput::new( population.get(), successes.get(), sample.get(), sample_successes.get(), ) + .map(hyper_geometric) .unwrap_or_default() }; view! { @@ -119,11 +120,12 @@ pub fn BinomCalculator() -> impl IntoView { let (trials_number, set_trials_number) = signal(0u8); let (successes_number, set_successes_number) = signal(0u8); let result = move || { - binomial( + BinomialInput::new( success_probability.get(), trials_number.get(), successes_number.get(), ) + .map(binomial) .unwrap_or_default() }; view! { @@ -135,6 +137,7 @@ pub fn BinomCalculator() -> impl IntoView { type="number" min=0 max=1 + prop:step=0.1 prop:value=success_probability on:input:target=move |ev| { set_success_probability.set(ev.target().value().parse().unwrap_or_default())