diff --git a/Cargo.lock b/Cargo.lock index c98132d..2e0f41e 100644 --- a/Cargo.lock +++ b/Cargo.lock @@ -20,6 +20,7 @@ dependencies = [ "csv", "env_logger", "fast-math", + "half", "lazy_static", "log", "plotters", @@ -200,6 +201,12 @@ version = "0.8.21" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "d0a5c400df2834b80a4c3327b3aad3a4c4cd4de0629063962b03235697506a28" +[[package]] +name = "crunchy" +version = "0.2.4" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "460fbee9c2c2f33933d720630a6a0bac33ba7053db5344fac858d4b8952d77d5" + [[package]] name = "crypto-common" version = "0.1.7" @@ -300,6 +307,18 @@ dependencies = [ "wasi", ] +[[package]] +name = "half" +version = "2.7.1" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "6ea2d84b969582b4b1864a92dc5d27cd2b77b622a8d79306834f1be5ba20d84b" +dependencies = [ + "cfg-if", + "crunchy", + "serde", + "zerocopy", +] + [[package]] name = "iana-time-zone" version = "0.1.63" @@ -804,18 +823,18 @@ dependencies = [ [[package]] name = "zerocopy" -version = "0.8.24" +version = "0.8.52" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "2586fea28e186957ef732a5f8b3be2da217d65c5969d4b1e17f973ebbe876879" +checksum = "ce1022995ff5ff5d841ad7d994facc23098cd40152f2c1d11cd607c6f530653f" dependencies = [ "zerocopy-derive", ] [[package]] name = "zerocopy-derive" -version = "0.8.24" +version = "0.8.52" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "a996a8f63c5c4448cd959ac1bab0aaa3306ccfd060472f85943ee0750f0169be" +checksum = "1ae7f38b72ec2a254e2b87ef277cf2cd4fb97cbebf944faa6f33354da0867930" dependencies = [ "proc-macro2", "quote", diff --git a/Cargo.toml b/Cargo.toml index 8e96ab6..31a825d 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -33,3 +33,4 @@ log = "0.4" env_logger = "0.11" fast-math = "0.1.1" sha2 = "0.10" +half = { version = "2.4", features = ["serde"] } diff --git a/src/rules/mod.rs b/src/rules/mod.rs index 9262117..42f1d07 100644 --- a/src/rules/mod.rs +++ b/src/rules/mod.rs @@ -91,9 +91,9 @@ use crate::config::{ RUN_PATHWAY_MICROBIOME_DISRUPTION_MULTIPLIER_KEY, RUN_PATHWAY_REVERSION_RATE_MULTIPLIER_KEY, }; use crate::simulation::population::{ - self, CarriageCompartment, HospitalStatus, ImmunodeficiencyType, Individual, - InfectionResolutionType, Region, ResistanceMechanism, BACTERIA_LIST, DRUG_CLASS_LOOKUP, - DRUG_SHORT_NAMES, INFECTION_EPS, MICROBIOME_MAJORITY_THRESHOLD, + self, load_float, store_float, CarriageCompartment, HospitalStatus, ImmunodeficiencyType, + Individual, InfectionResolutionType, Region, ResistanceMechanism, BACTERIA_LIST, + DRUG_CLASS_LOOKUP, DRUG_SHORT_NAMES, INFECTION_EPS, MICROBIOME_MAJORITY_THRESHOLD, }; use rand::Rng; @@ -311,7 +311,7 @@ fn has_serious_resistance_test_positive(individual: &Individual) -> bool { for &drug_name in serious_resistance_marker_drugs(bacteria_name) { if let Some(&drug_idx) = DRUG_INDEX_BY_NAME.get(drug_name) { - if individual.resistances[b_idx][drug_idx].test_r > INFECTION_EPS { + if load_float(individual.resistances[b_idx][drug_idx].test_r) > INFECTION_EPS { return true; } } @@ -443,12 +443,12 @@ fn propagate_mechanism_resistance( if raise_only { // Acquisition mode only raises resistance: a cache-sampled resistant profile // should not be weakened by a partial set of currently visible mechanisms. - if new_any_r > resistance_data.any_r { - resistance_data.any_r = new_any_r; + if new_any_r > load_float(resistance_data.any_r) { + resistance_data.any_r = store_float(new_any_r); } } else { // Reset mode (reversion) — set any_r to exact mechanism-derived level - resistance_data.any_r = new_any_r; + resistance_data.any_r = store_float(new_any_r); } // Derive microbiome_r from mechanism_microbiome state @@ -457,11 +457,11 @@ fn propagate_mechanism_resistance( .min(max_resistance_level) .max(0.0); if raise_only { - if new_microbiome_r > resistance_data.microbiome_r { - resistance_data.microbiome_r = new_microbiome_r; + if new_microbiome_r > load_float(resistance_data.microbiome_r) { + resistance_data.microbiome_r = store_float(new_microbiome_r); } } else { - resistance_data.microbiome_r = new_microbiome_r; + resistance_data.microbiome_r = store_float(new_microbiome_r); } } } @@ -2740,7 +2740,7 @@ pub(crate) fn apply_rules( // TMP-SMX is valid for PCP prophylaxis but shouldn't dominate the pool. // ^^^^ "trim_sulf" => 0.8, - "azithromycin" => 3.5, // 4.5 + "azithromycin" => 3.5, // 4.5 "ciprofloxacin" => 2.0, "levofloxacin" => 1.5, _ => 0.0, @@ -2764,7 +2764,7 @@ pub(crate) fn apply_rules( if (time_step as i32) < initiated_day + resistance_test_result_delay_days { continue; } - if individual.resistances[b_idx][drug_idx].test_r <= 0.0 { + if load_float(individual.resistances[b_idx][drug_idx].test_r) <= 0.0 { continue; } resistance_detected = true; @@ -3670,7 +3670,7 @@ pub(crate) fn apply_rules( // if we don't have a specific sensitivity test result confirming susceptibility. let sensitivity_result_available = identified_bacteria .iter() - .any(|&b_idx| individual.resistances[b_idx][drug_idx].test_r > 0.0); + .any(|&b_idx| load_float(individual.resistances[b_idx][drug_idx].test_r) > 0.0); if !sensitivity_result_available { let mut regional_resistance_penalty = 1.0_f64; @@ -4758,12 +4758,14 @@ pub(crate) fn apply_rules( let mut total_log_odds = base_log_odds; // (1) Bacteria level effect - higher bacteria level decreases recovery probability - let bacteria_level_coefficient = store.globals.sepsis_recovery_log_odds_bacteria_level; + let bacteria_level_coefficient = + store.globals.sepsis_recovery_log_odds_bacteria_level; total_log_odds += individual.level[b_idx] * bacteria_level_coefficient; // (2) Hospital status effect - being in hospital increases recovery probability if individual.hospital_status.is_hospitalized() { - let hospital_coefficient = store.globals.sepsis_recovery_log_odds_in_hospital; + let hospital_coefficient = + store.globals.sepsis_recovery_log_odds_in_hospital; total_log_odds += hospital_coefficient; } @@ -4788,7 +4790,9 @@ pub(crate) fn apply_rules( } // (5) Region-specific effect (healthcare quality and ICU availability) - total_log_odds += store.region.sepsis_recovery_log_odds(individual.region_living); + total_log_odds += store + .region + .sepsis_recovery_log_odds(individual.region_living); // Convert log odds to probability using logistic function let recovery_probability = 1.0 / (1.0 + (-total_log_odds).fast_exp()); @@ -5085,7 +5089,7 @@ pub(crate) fn apply_rules( individual.microbiome_acquired_on_drug_today[b_idx] = false; individual.microbiome_cleared_today[b_idx] = false; for resistance_data in individual.resistances[b_idx].iter_mut() { - resistance_data.microbiome_r = 0.0; + resistance_data.microbiome_r = store_float(0.0); } individual.clear_microbiome_mechanisms(b_idx); } @@ -5145,7 +5149,8 @@ pub(crate) fn apply_rules( let normalized_micro_r = if max_resistance_level <= f64::EPSILON { 1.0 } else { - (resistance_data.microbiome_r / max_resistance_level).clamp(0.0, 1.0) + (load_float(resistance_data.microbiome_r) / max_resistance_level) + .clamp(0.0, 1.0) }; let base_potency = param_cache.potency(b_idx, d_idx); let effective_activity = @@ -5330,7 +5335,7 @@ pub(crate) fn apply_rules( } else if !individual.presence_microbiome[b_idx] { // No microbiome presence — ensure microbiome_r is zero for d_idx in 0..DRUG_SHORT_NAMES.len() { - individual.resistances[b_idx][d_idx].microbiome_r = 0.0; + individual.resistances[b_idx][d_idx].microbiome_r = store_float(0.0); } } @@ -5351,7 +5356,8 @@ pub(crate) fn apply_rules( let base_potency = param_cache.potency(b_idx, drug_idx); let drug_current_level = individual.cur_level_drug[drug_idx]; let max_resistance_level = store.globals.max_resistance_level; - let resistance_level = individual.resistances[b_idx][drug_idx].any_r; + let resistance_level = + load_float(individual.resistances[b_idx][drug_idx].any_r); let normalized_any_r = resistance_level / max_resistance_level; let effective_activity = base_potency * drug_current_level * (1.0 - normalized_any_r); @@ -5723,7 +5729,7 @@ pub(crate) fn apply_rules( // Provenance bookkeeping disabled for memory-saving runs. if crate::simulation::population::TRACK_RESISTANCE_ACQUISITION_PROVENANCE { for d_idx in 0..DRUG_SHORT_NAMES.len() { - if individual.resistances[b_idx][d_idx].any_r > 0.0 + if load_float(individual.resistances[b_idx][d_idx].any_r) > 0.0 && individual.how_resistance_acquired[b_idx][d_idx].is_none() { individual.how_resistance_acquired[b_idx][d_idx] = Some( @@ -5941,8 +5947,11 @@ pub(crate) fn apply_rules( &mut individual.resistances[bacteria_full_idx][drug_index]; // Clamp any_r to valid range - resistance_data.any_r = - resistance_data.any_r.min(max_resistance_level).max(0.0); + resistance_data.any_r = store_float( + load_float(resistance_data.any_r) + .min(max_resistance_level) + .max(0.0), + ); if drug_current_level > 0.0 { // Fetch potency from cached lookup @@ -5950,7 +5959,8 @@ pub(crate) fn apply_rules( // any_r is updated when mechanism state changes; this loop only // translates the current resistance state into per-drug activity. - let normalized_any_r = resistance_data.any_r / max_resistance_level; + let normalized_any_r = + load_float(resistance_data.any_r) / max_resistance_level; // Apply syndrome-specific drug penetration factor // This accounts for pharmacokinetic differences at different infection sites @@ -5964,10 +5974,11 @@ pub(crate) fn apply_rules( // Effective drug level at infection site let effective_drug_level = drug_current_level * penetration_factor; - resistance_data.activity_r = - base_potency * effective_drug_level * (1.0 - normalized_any_r); + resistance_data.activity_r = store_float( + base_potency * effective_drug_level * (1.0 - normalized_any_r), + ); } else { - resistance_data.activity_r = 0.0; + resistance_data.activity_r = store_float(0.0); } } } @@ -6047,13 +6058,14 @@ pub(crate) fn apply_rules( if test_initiated_day != -1 && (time_step as i32) >= (test_initiated_day + resistance_test_result_delay_days) { - let test_r_already_set = - individual.resistances[b_idx].iter().any(|r| r.test_r > 0.0); + let test_r_already_set = individual.resistances[b_idx] + .iter() + .any(|r| load_float(r.test_r) > 0.0); if !test_r_already_set { for d_idx in 0..DRUG_SHORT_NAMES.len() { // Use any_r to model standard clinical microbiologic phenotypic testing. // (majority_r removed — any_r is the primary resistance measure) - let major_r = individual.resistances[b_idx][d_idx].any_r; + let major_r = load_float(individual.resistances[b_idx][d_idx].any_r); let error = rng.gen_bool(test_r_error_prob); let test_r = if error { if major_r < INFECTION_EPS { @@ -6064,14 +6076,14 @@ pub(crate) fn apply_rules( } else { major_r }; - individual.resistances[b_idx][d_idx].test_r = test_r; + individual.resistances[b_idx][d_idx].test_r = store_float(test_r); } } } } else { // Reset resistance test results if bacterial identification test is negative for d_idx in 0..DRUG_SHORT_NAMES.len() { - individual.resistances[b_idx][d_idx].test_r = 0.0; + individual.resistances[b_idx][d_idx].test_r = store_float(0.0); } } @@ -6210,7 +6222,8 @@ pub(crate) fn apply_rules( || individual.microbiome_mechanism_mask(b_idx) != 0; if !has_active_mechanism && individual.resistances[b_idx].iter().any(|resistance| { - resistance.any_r > 0.0 || resistance.microbiome_r > 0.0 + load_float(resistance.any_r) > 0.0 + || load_float(resistance.microbiome_r) > 0.0 }) { let reversion_rate = store @@ -6223,10 +6236,10 @@ pub(crate) fn apply_rules( for drug_index in 0..DRUG_SHORT_NAMES.len() { let resistance_data = &mut individual.resistances[b_idx][drug_index]; - resistance_data.microbiome_r = 0.0; - resistance_data.test_r = 0.0; - resistance_data.activity_r = 0.0; - resistance_data.any_r = 0.0; + resistance_data.microbiome_r = store_float(0.0); + resistance_data.test_r = store_float(0.0); + resistance_data.activity_r = store_float(0.0); + resistance_data.any_r = store_float(0.0); // Provenance bookkeeping disabled for memory-saving runs. if crate::simulation::population::TRACK_RESISTANCE_ACQUISITION_PROVENANCE { individual.how_resistance_acquired[b_idx][drug_index] = None; @@ -6240,7 +6253,7 @@ pub(crate) fn apply_rules( for (drug_idx, _drug_name) in DRUG_SHORT_NAMES.iter().enumerate() { if individual.cur_level_drug[drug_idx] > 0.0 { let resistance_data = &individual.resistances[b_idx][drug_idx]; - total_reduction_due_to_antibiotic += resistance_data.activity_r; + total_reduction_due_to_antibiotic += load_float(resistance_data.activity_r); } } @@ -6326,7 +6339,7 @@ pub(crate) fn apply_rules( ) { let activity_values: Vec = individual.resistances[b_idx] .iter() - .map(|resistance| resistance.activity_r) + .map(|resistance| load_float(resistance.activity_r)) .collect(); crate::simulation::journey_logger::cache_pre_clearance_activity( individual.id, @@ -6364,7 +6377,7 @@ pub(crate) fn apply_rules( if individual.resistances[b_idx] .iter() - .any(|resistance| resistance.any_r > 0.0) + .any(|resistance| load_float(resistance.any_r) > 0.0) { let category = individual .microbiome_resistance_level(b_idx, MICROBIOME_MAJORITY_THRESHOLD); @@ -6390,8 +6403,8 @@ pub(crate) fn apply_rules( // Clear infection data after tracking resolution for drug_idx_clear in 0..DRUG_SHORT_NAMES.len() { let resistance_data = &mut individual.resistances[b_idx][drug_idx_clear]; - resistance_data.any_r = 0.0; - resistance_data.activity_r = 0.0; + resistance_data.any_r = store_float(0.0); + resistance_data.activity_r = store_float(0.0); // Provenance bookkeeping disabled for memory-saving runs. if crate::simulation::population::TRACK_RESISTANCE_ACQUISITION_PROVENANCE { individual.how_resistance_acquired[b_idx][drug_idx_clear] = None; @@ -6505,8 +6518,9 @@ pub(crate) fn apply_rules( compartment_masks[b_idx] = mask; infection_presence[b_idx] = has_infection; - let has_any_resistance = - individual.resistances[b_idx].iter().any(|r| r.any_r > 0.0); + let has_any_resistance = individual.resistances[b_idx] + .iter() + .any(|r| load_float(r.any_r) > 0.0); if has_any_resistance { potential_donors.push(b_idx); @@ -6638,7 +6652,8 @@ pub(crate) fn apply_rules( // Provenance bookkeeping disabled for memory-saving runs. if crate::simulation::population::TRACK_RESISTANCE_ACQUISITION_PROVENANCE { for drug_idx in 0..DRUG_SHORT_NAMES.len() { - if individual.resistances[recipient_idx][drug_idx].any_r > 0.0 + if load_float(individual.resistances[recipient_idx][drug_idx].any_r) + > 0.0 && individual.how_resistance_acquired[recipient_idx][drug_idx] .is_none() { @@ -6705,8 +6720,9 @@ fn apply_cross_resistance( if let Some(resistance_data) = individual.resistances.get(b_idx).and_then(|r| r.get(d_idx)) { - if resistance_data.any_r > max_any_r { - max_any_r = resistance_data.any_r; + let any_r = load_float(resistance_data.any_r); + if any_r > max_any_r { + max_any_r = any_r; } } } @@ -6719,7 +6735,7 @@ fn apply_cross_resistance( .get_mut(b_idx) .and_then(|r| r.get_mut(d_idx)) { - resistance_data.any_r = max_any_r; + resistance_data.any_r = store_float(max_any_r); } } } diff --git a/src/simulation/journey_logger.rs b/src/simulation/journey_logger.rs index 9ad2836..beb1e75 100644 --- a/src/simulation/journey_logger.rs +++ b/src/simulation/journey_logger.rs @@ -6,7 +6,7 @@ // resistance, and clearance trajectories for a single journey. use crate::simulation::population::{ - Individual, ResistanceMechanism, BACTERIA_LIST, DRUG_SHORT_NAMES, INFECTION_EPS, + load_float, Individual, ResistanceMechanism, BACTERIA_LIST, DRUG_SHORT_NAMES, INFECTION_EPS, MICROBIOME_MAJORITY_THRESHOLD, }; use rand::rngs::SmallRng; @@ -652,7 +652,8 @@ impl JourneyLogger { .iter() .enumerate() .filter(|(idx, _)| { - let resistance_value = individual.resistances[primary_bacteria_idx][*idx].any_r; + let resistance_value = + load_float(individual.resistances[primary_bacteria_idx][*idx].any_r); let drug_active = individual.cur_use_drug[*idx] || individual.cur_level_drug[*idx] > INFECTION_EPS; resistance_value > 0.0 || drug_active @@ -660,7 +661,7 @@ impl JourneyLogger { .map(|(idx, &drug_name)| { ( drug_name.to_string(), - individual.resistances[primary_bacteria_idx][idx].any_r, + load_float(individual.resistances[primary_bacteria_idx][idx].any_r), ) }) .collect(); @@ -684,7 +685,7 @@ impl JourneyLogger { .map(|(idx, &drug_name)| { ( drug_name.to_string(), - individual.resistances[primary_bacteria_idx][idx].activity_r, + load_float(individual.resistances[primary_bacteria_idx][idx].activity_r), ) }) .collect(); @@ -722,7 +723,8 @@ impl JourneyLogger { let mut resistances_microbiome_r: Vec<(String, f64)> = Vec::new(); for (idx, &drug_name) in DRUG_SHORT_NAMES.iter().enumerate() { - let microbiome_value = individual.resistances[primary_bacteria_idx][idx].microbiome_r; + let microbiome_value = + load_float(individual.resistances[primary_bacteria_idx][idx].microbiome_r); if primary_microbiome_present || microbiome_value > 0.0 { resistances_microbiome_r.push((drug_name.to_string(), microbiome_value)); } @@ -741,7 +743,7 @@ impl JourneyLogger { let mut majority_max = 0.0; for resistance in resistances { - let value = resistance.microbiome_r; + let value = load_float(resistance.microbiome_r); if value >= MICROBIOME_MAJORITY_THRESHOLD { if value > majority_max { majority_max = value; @@ -1226,7 +1228,8 @@ impl JourneyLogger { .find(|(name, _)| name == drug_name) .map(|(_, value)| *value) .unwrap_or(0.0); - let current_any = individual.resistances[primary_bacteria_idx][drug_idx].any_r; + let current_any = + load_float(individual.resistances[primary_bacteria_idx][drug_idx].any_r); let any_increased = current_any > prev_any + RESISTANCE_EPSILON; diff --git a/src/simulation/population.rs b/src/simulation/population.rs index 9327463..b5dd03f 100644 --- a/src/simulation/population.rs +++ b/src/simulation/population.rs @@ -65,6 +65,74 @@ use std::fmt; // Minimum infection/drug level threshold; values below this are treated as cleared to avoid floating-point noise. pub const INFECTION_EPS: f64 = 0.001; +pub trait StoredFloatRepr: Copy { + fn store(value: f64) -> Self; + fn load(self) -> f64; +} + +impl StoredFloatRepr for f64 { + #[inline] + fn store(value: f64) -> Self { + value + } + + #[inline] + fn load(self) -> f64 { + self + } +} + +impl StoredFloatRepr for half::f16 { + #[inline] + fn store(value: f64) -> Self { + half::f16::from_f64(value) + } + + #[inline] + fn load(self) -> f64 { + self.to_f64() + } +} + +#[derive(Debug, Copy, Clone, Serialize, Deserialize)] +pub struct BoundedResistanceU16(u16); + +impl StoredFloatRepr for BoundedResistanceU16 { + #[inline] + fn store(value: f64) -> Self { + let bounded = if value.is_finite() { + value.clamp(0.0, 1.0) + } else { + 0.0 + }; + BoundedResistanceU16((bounded * u16::MAX as f64).round() as u16) + } + + #[inline] + fn load(self) -> f64 { + self.0 as f64 / u16::MAX as f64 + } +} + +impl fmt::Display for BoundedResistanceU16 { + fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result { + self.load().fmt(f) + } +} + +pub type StoredBoundedResistanceFloat = BoundedResistanceU16; +pub type StoredActivityResistanceFloat = half::f16; + +#[inline] +pub fn store_float(value: f64) -> T { + T::store(value) +} + +#[inline] +pub fn load_float(value: T) -> f64 { + value.load() +} + /// Specific resistance mechanisms that can be present in bacteria /// These mechanism bitmasks are the source of infection, majority-strain, and /// microbiome resistance projections. @@ -1254,20 +1322,20 @@ impl Region { pub struct Resistance { /// Resistance level in colonizing (carriage) bacteria. Range: 0.0-1.0. /// Derived from mechanism_microbiome via the multiplicative product formula. - pub microbiome_r: f64, + pub microbiome_r: StoredBoundedResistanceFloat, /// Resistance as would be detected by laboratory testing. Range: 0.0-1.0. /// May differ from actual resistance due to test sensitivity/specificity. - pub test_r: f64, + pub test_r: StoredBoundedResistanceFloat, /// Effective resistance for drug activity calculations. Range: 0.0-1.0. /// Takes into account mechanism-specific effects on drug binding/activity. - pub activity_r: f64, + pub activity_r: StoredActivityResistanceFloat, /// Primary resistance level - resistance present in ANY infected bacteria. Range: 0.0-1.0. /// Derived from mechanism_any via the multiplicative product formula: /// any_r = 1 - product(1 - enhancement_multiplier) over all present mechanisms. - pub any_r: f64, + pub any_r: StoredBoundedResistanceFloat, } pub const MICROBIOME_RESISTANCE_LEVEL_COUNT: usize = 4; @@ -1787,10 +1855,10 @@ impl Individual { let mut drug_resistances = Vec::with_capacity(num_drugs); for _ in 0..num_drugs { drug_resistances.push(Resistance { - microbiome_r: 0.0, - test_r: 0.0, - activity_r: 0.0, - any_r: 0.0, + microbiome_r: store_float(0.0), + test_r: store_float(0.0), + activity_r: store_float(0.0), + any_r: store_float(0.0), }); } resistances.push(drug_resistances); @@ -1935,9 +2003,10 @@ impl Individual { if let Some(resistances) = self.resistances.get(bacteria_idx) { for resistance in resistances { - if resistance.microbiome_r > 0.0 { + let microbiome_r = load_float(resistance.microbiome_r); + if microbiome_r > 0.0 { has_resistance = true; - if resistance.microbiome_r >= threshold { + if microbiome_r >= threshold { has_majority = true; break; } diff --git a/src/simulation/simulation.rs b/src/simulation/simulation.rs index 76b3695..1b1c95a 100644 --- a/src/simulation/simulation.rs +++ b/src/simulation/simulation.rs @@ -16,8 +16,8 @@ use crate::observability; use crate::rules::apply_rules; use crate::simulation::journey_logger::JourneyLogger; use crate::simulation::population::{ - MicrobiomeResistanceLevel, Population, Region, ResistanceMechanism, BACTERIA_LIST, - DRUG_SHORT_NAMES, INFECTION_EPS, MICROBIOME_MAJORITY_THRESHOLD, + load_float, store_float, MicrobiomeResistanceLevel, Population, Region, ResistanceMechanism, + BACTERIA_LIST, DRUG_SHORT_NAMES, INFECTION_EPS, MICROBIOME_MAJORITY_THRESHOLD, MICROBIOME_RESISTANCE_LEVEL_COUNT, }; use crate::simulation::rng::{ @@ -1286,10 +1286,10 @@ impl IndividualLogger { let mut any_r = Vec::new(); for bact in &ind.resistances { for res in bact { - microbiome_r.push(res.microbiome_r); - test_r.push(res.test_r); - activity_r.push(res.activity_r); - any_r.push(res.any_r); + microbiome_r.push(load_float(res.microbiome_r)); + test_r.push(load_float(res.test_r)); + activity_r.push(load_float(res.activity_r)); + any_r.push(load_float(res.any_r)); } } @@ -1317,7 +1317,7 @@ impl IndividualLogger { if ind.level[b_idx] > 0.0 && ind.cur_use_drug.iter().any(|&on_drug| on_drug) { for d_idx in 0..DRUG_SHORT_NAMES.len() { if ind.cur_use_drug[d_idx] { - result = ind.resistances[b_idx][d_idx].activity_r; + result = load_float(ind.resistances[b_idx][d_idx].activity_r); break; } } @@ -3359,22 +3359,22 @@ impl Simulation { } for d_idx in 0..num_drugs { let resistance_data = &individual.resistances[b_idx][d_idx]; + let any_r = load_float(resistance_data.any_r); if need_full_summary { let threshold = mic_lt2_thresholds[base + d_idx]; - if resistance_data.any_r < threshold { + if any_r < threshold { lt.mic_lt2_counts[base + d_idx] += 1; } - lt.any_r_sum_by_bacteria_drug[base + d_idx] += resistance_data.any_r; + lt.any_r_sum_by_bacteria_drug[base + d_idx] += any_r; let potency = potency_matrix[base + d_idx]; let mic = if potency <= 1e-9 { 1e12 } else { - let susceptible_fraction = - (1.0 - resistance_data.any_r).clamp(1e-6, 1.0); + let susceptible_fraction = (1.0 - any_r).clamp(1e-6, 1.0); 1.0 / (susceptible_fraction * potency) }; lt.mic_sum_by_bacteria_drug[base + d_idx] += mic; - if resistance_data.any_r > 0.0 { + if any_r > 0.0 { lt.infected_with_any_r_positive_by_bacteria_drug[base + d_idx] += 1; if !lt.infected_with_any_r_positive_hospital_by_bacteria_drug.is_empty() { if individual.hospital_status.is_hospitalized() { @@ -3385,11 +3385,11 @@ impl Simulation { } } if individual.infection_hospital_acquired[b_idx] { - lt.any_r_sum_by_bacteria_drug_hospital[base + d_idx] += resistance_data.any_r; + lt.any_r_sum_by_bacteria_drug_hospital[base + d_idx] += any_r; } } if let Some(region_idx) = effective_region_idx_for_any_r { - lt.any_r_sum_by_region[region_idx] += resistance_data.any_r; + lt.any_r_sum_by_region[region_idx] += any_r; } } @@ -3429,7 +3429,7 @@ impl Simulation { if has_any_microbiome && need_full_summary { for d_idx in 0..num_drugs { let resistance_data = &individual.resistances[b_idx][d_idx]; - if resistance_data.microbiome_r > 0.0 { + if load_float(resistance_data.microbiome_r) > 0.0 { let idx = b_idx * num_drugs + d_idx; lt.microbiome_r_positive_by_bacteria_drug[idx] += 1; } @@ -3743,7 +3743,7 @@ impl Simulation { } let has_resistant_microbiome = individual.resistances[b_idx] .iter() - .any(|resistance| resistance.microbiome_r > 0.0); + .any(|resistance| load_float(resistance.microbiome_r) > 0.0); if has_resistant_microbiome { lt.presence_microbiome_resistant_by_bacteria[b_idx] += 1; } @@ -3930,7 +3930,9 @@ impl Simulation { lt.newly_infected_hospital_by_bacteria_region[b_idx * 6 + cur_region_idx] += 1; } // Resistance-at-infection tracking - let has_any_r = individual.resistances[b_idx].iter().any(|rd| rd.any_r > 0.0); + let has_any_r = individual.resistances[b_idx] + .iter() + .any(|rd| load_float(rd.any_r) > 0.0); if has_any_r && !lt.newly_infected_any_r_hospital_by_bacteria.is_empty() { if individual.infection_hospital_acquired[b_idx] { lt.newly_infected_any_r_hospital_by_bacteria[b_idx] += 1; @@ -3944,25 +3946,26 @@ impl Simulation { // Full iteration for stats that need all drugs. for d_idx in 0..num_drugs { let resistance_data = &individual.resistances[b_idx][d_idx]; + let any_r = load_float(resistance_data.any_r); // Only sum activity_r if individual is currently on this drug if individual.cur_use_drug[d_idx] { - activity_r_sum += resistance_data.activity_r; + activity_r_sum += load_float(resistance_data.activity_r); let base_potency = self.param_cache.potency(b_idx, d_idx); let syndrome_id = individual.infectious_syndrome[b_idx] as usize; let penetration_factor = config::parameter_store().syndrome.drug_penetration(syndrome_id, d_idx); let effective_drug_level = individual.cur_level_drug[d_idx] * penetration_factor; - let normalized_any_r = resistance_data.any_r - / config::parameter_store().globals.max_resistance_level; + let normalized_any_r = + any_r / config::parameter_store().globals.max_resistance_level; max_possible_activity_r_sum += base_potency * effective_drug_level; lt.activity_r_pure_sum_by_bacteria[b_idx] += base_potency * (1.0 - normalized_any_r); lt.max_possible_activity_r_pure_sum_by_bacteria[b_idx] += base_potency; } - if need_full_summary && resistance_data.any_r > 0.0 { + if need_full_summary && any_r > 0.0 { lt.resistance_by_bacteria_drug[base + d_idx] += 1; } - if resistance_data.any_r > 0.0 { + if any_r > 0.0 { infection_any_r_positive = true; individual_has_any_r_positive = true; @@ -4845,10 +4848,10 @@ impl Simulation { for b_idx in 0..num_bacteria { for d_idx in 0..num_drugs { let resistance = &mut individual.resistances[b_idx][d_idx]; - resistance.any_r = 0.0; - resistance.activity_r = 0.0; - resistance.microbiome_r = 0.0; - resistance.test_r = 0.0; + resistance.any_r = store_float(0.0); + resistance.activity_r = store_float(0.0); + resistance.microbiome_r = store_float(0.0); + resistance.test_r = store_float(0.0); // Provenance bookkeeping is disabled for memory-saving calibration // runs, so there may be no dense provenance matrix to clear here. if crate::simulation::population::TRACK_RESISTANCE_ACQUISITION_PROVENANCE { diff --git a/tests/csv_invariants.rs b/tests/csv_invariants.rs index 86152e6..981f998 100644 --- a/tests/csv_invariants.rs +++ b/tests/csv_invariants.rs @@ -1,5 +1,5 @@ use amr_project::simulation::journey_logger::JourneyLogger; -use amr_project::simulation::population::Individual; +use amr_project::simulation::population::{store_float, Individual}; use amr_project::simulation::simulation::{CalibrationMode, Simulation}; use rand::rngs::SmallRng; use rand::SeedableRng; @@ -99,9 +99,9 @@ fn journey_csv_rows_match_header_width_for_forced_infection() { individual.cur_level_drug[drug_idx] = 5.0; individual.date_drug_initiated[drug_idx] = 1; individual.days_on_current_treatment[bacteria_idx] = 1; - individual.resistances[bacteria_idx][drug_idx].any_r = 0.25; - individual.resistances[bacteria_idx][drug_idx].activity_r = 0.75; - individual.resistances[bacteria_idx][drug_idx].microbiome_r = 0.50; + individual.resistances[bacteria_idx][drug_idx].any_r = store_float(0.25); + individual.resistances[bacteria_idx][drug_idx].activity_r = store_float(0.75); + individual.resistances[bacteria_idx][drug_idx].microbiome_r = store_float(0.50); let mut logger = JourneyLogger::new(Some(123_456_789)); logger