Skip to content

Commit fdadd61

Browse files
sirnosirno
committed
Big commit, using new Symbol encoding
1 parent 574438f commit fdadd61

24 files changed

+566
-532
lines changed

‎CHANGELOG.md

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -10,6 +10,7 @@
1010
- Huge performance improvements when collecting all mutations
1111
- Tree will now merge internal nodes when possible
1212
- Methods that access data have to be annotated with `#[require_deferred_drop]`
13+
- Introduced generic encoding for symbolic values `virolution::encoding::Symbol`
1314

1415
## 0.3.0 --- Mixed Performance (May 06, 2024)
1516

‎Cargo.lock

Lines changed: 12 additions & 0 deletions
Some generated files are not rendered by default. Learn more about customizing how changed files appear on GitHub.

‎Cargo.toml

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -30,6 +30,7 @@ clap = { version = "4.3.11", features = ["derive"] }
3030
csv = "1.2.2"
3131
derivative = "2.2.0"
3232
derive_more = "0.99.17"
33+
derive-where = "1.2.7"
3334
evalexpr = "8.2.0"
3435
indicatif = "0.17.5"
3536
itertools = "0.10.5"

‎examples/compartments.rs

Lines changed: 5 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -8,14 +8,15 @@ use virolution::core::fitness::utility::UtilityFunction;
88
use virolution::core::fitness::FitnessProvider;
99
use virolution::core::haplotype::*;
1010
use virolution::core::population::Population;
11+
use virolution::encoding::Nucleotide as Nt;
1112
use virolution::simulation::*;
1213

1314
use virolution::population;
1415

1516
fn main() {
1617
let plan_path = PathBuf::from_iter([env!("CARGO_MANIFEST_DIR"), "data/plan.csv"]);
1718
let plan = Schedule::read(plan_path.to_str().unwrap()).expect("Failed to read plan");
18-
let sequence = vec![Some(0x00); 5386];
19+
let sequence = vec![Nt::A; 5386];
1920
let distribution = FitnessDistribution::Exponential(ExponentialParameters {
2021
weights: MutationCategoryWeights {
2122
beneficial: 0.29,
@@ -28,7 +29,7 @@ fn main() {
2829
});
2930
let fitness_model = FitnessModel::new(distribution.clone(), UtilityFunction::Linear);
3031

31-
let fitness_provider = FitnessProvider::from_model(0, &sequence, 4, &fitness_model)
32+
let fitness_provider = FitnessProvider::from_model(0, &sequence, &fitness_model)
3233
.expect("Failed to create fitness table");
3334

3435
let wt = Wildtype::new(sequence);
@@ -50,7 +51,7 @@ fn main() {
5051
};
5152

5253
let n_compartments = 3;
53-
let mut compartment_simulations: Vec<BasicSimulation> = (0..n_compartments)
54+
let mut compartment_simulations: Vec<BasicSimulation<Nt>> = (0..n_compartments)
5455
.map(|_| {
5556
let population = population![wt.clone(); 1_000_000];
5657
let fitness_providers =
@@ -78,8 +79,7 @@ fn main() {
7879

7980
let transfers = plan.get_transfer_matrix(gen);
8081

81-
let populations: Vec<Population> = (0..n_compartments)
82-
.into_iter()
82+
let populations: Vec<Population<Nt>> = (0..n_compartments)
8383
.map(|target| {
8484
Population::from_iter((0..n_compartments).map(|origin| {
8585
compartment_simulations[origin]

‎examples/haplotype.rs

Lines changed: 10 additions & 15 deletions
Original file line numberDiff line numberDiff line change
@@ -1,13 +1,14 @@
11
extern crate virolution;
22

33
use virolution::core::haplotype::*;
4+
use virolution::encoding::Nucleotide as Nt;
45

56
fn main() {
6-
let bytes = vec![Some(0x00); 4];
7+
let bytes = vec![Nt::A; 4];
78
let wt = Wildtype::new(bytes);
8-
let ht = wt.create_descendant(vec![2], vec![Some(0x01)], 0);
9-
let ht2 = ht.create_descendant(vec![1], vec![Some(0x02)], 0);
10-
let ht3 = ht2.create_descendant(vec![2], vec![Some(0x03)], 0);
9+
let ht = wt.create_descendant(vec![2], vec![Nt::T], 0);
10+
let ht2 = ht.create_descendant(vec![1], vec![Nt::C], 0);
11+
let ht3 = ht2.create_descendant(vec![2], vec![Nt::G], 0);
1112
let ht4 = Haplotype::create_recombinant(&ht, &ht3, 0, 2, 0);
1213
println!("---debug---");
1314
println!("wt: {:?}", *wt);
@@ -18,29 +19,23 @@ fn main() {
1819
println!("---get_base---");
1920
println!(
2021
"{:?}",
21-
(0..4).map(|idx| wt.get_base(&idx)).collect::<Vec<Symbol>>()
22+
(0..4).map(|idx| wt.get_base(&idx)).collect::<Vec<Nt>>()
2223
);
2324
println!(
2425
"{:?}",
25-
(0..4).map(|idx| ht.get_base(&idx)).collect::<Vec<Symbol>>()
26+
(0..4).map(|idx| ht.get_base(&idx)).collect::<Vec<Nt>>()
2627
);
2728
println!(
2829
"{:?}",
29-
(0..4)
30-
.map(|idx| ht2.get_base(&idx))
31-
.collect::<Vec<Symbol>>()
30+
(0..4).map(|idx| ht2.get_base(&idx)).collect::<Vec<Nt>>()
3231
);
3332
println!(
3433
"{:?}",
35-
(0..4)
36-
.map(|idx| ht3.get_base(&idx))
37-
.collect::<Vec<Symbol>>()
34+
(0..4).map(|idx| ht3.get_base(&idx)).collect::<Vec<Nt>>()
3835
);
3936
println!(
4037
"{:?}",
41-
(0..4)
42-
.map(|idx| ht4.get_base(&idx))
43-
.collect::<Vec<Symbol>>()
38+
(0..4).map(|idx| ht4.get_base(&idx)).collect::<Vec<Nt>>()
4439
);
4540

4641
println!("---get_sequence---");

‎examples/simulation.rs

Lines changed: 7 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -7,12 +7,13 @@ use virolution::core::fitness::utility::UtilityFunction;
77
use virolution::core::fitness::FitnessProvider;
88
use virolution::core::haplotype::*;
99
use virolution::core::Population;
10+
use virolution::encoding::Nucleotide as Nt;
1011
use virolution::simulation::*;
1112

1213
use virolution::population;
1314

1415
fn main() {
15-
let sequence = vec![Some(0x00); 100];
16+
let sequence = vec![Nt::A; 100];
1617
let distribution = FitnessDistribution::Exponential(ExponentialParameters {
1718
weights: MutationCategoryWeights {
1819
beneficial: 0.29,
@@ -25,12 +26,12 @@ fn main() {
2526
});
2627
let fitness_model = FitnessModel::new(distribution.clone(), UtilityFunction::Linear);
2728

28-
let fitness_table = FitnessProvider::from_model(0, &sequence, 4, &fitness_model).unwrap();
29+
let fitness_table = FitnessProvider::from_model(0, &sequence, &fitness_model).unwrap();
2930

3031
let wt = Wildtype::new(sequence);
31-
let ht = wt.create_descendant(vec![2], vec![Some(0x01)], 0);
32-
let ht2 = ht.create_descendant(vec![1], vec![Some(0x02)], 0);
33-
let ht3 = ht2.create_descendant(vec![2], vec![Some(0x03)], 0);
32+
let ht = wt.create_descendant(vec![2], vec![Nt::T], 0);
33+
let ht2 = ht.create_descendant(vec![1], vec![Nt::C], 0);
34+
let ht3 = ht2.create_descendant(vec![2], vec![Nt::G], 0);
3435
let ht4 = Haplotype::create_recombinant(&ht, &ht3, 0, 2, 0);
3536

3637
println!("---fitnesses---");
@@ -40,7 +41,7 @@ fn main() {
4041
println!("ht3: {}", ht3.get_fitness(&fitness_table));
4142
println!("ht4: {}", ht4.get_fitness(&fitness_table));
4243

43-
let population: Population = population![wt.clone(); 10];
44+
let population: Population<Nt> = population![wt.clone(); 10];
4445
let simulation_settings = Parameters {
4546
mutation_rate: 1e-6,
4647
recombination_rate: 0.,

‎src/core/fitness/epistasis.rs

Lines changed: 19 additions & 18 deletions
Original file line numberDiff line numberDiff line change
@@ -1,13 +1,13 @@
11
use npyz::WriterBuilder;
22
use std::collections::HashMap;
33

4-
use crate::core::haplotype::Symbol;
4+
use crate::encoding::Symbol;
55
use crate::errors::VirolutionError;
66
use crate::references::HaplotypeRef;
77

88
use super::init::{FitnessDistribution, FitnessModel};
99

10-
type EpistasisTableKey = (usize, Symbol);
10+
type EpistasisTableKey<S> = (usize, S);
1111

1212
#[derive(Debug, npyz::Deserialize, npyz::Serialize, npyz::AutoSerialize)]
1313
pub struct EpiEntry {
@@ -19,11 +19,11 @@ pub struct EpiEntry {
1919
}
2020

2121
#[derive(Clone, Debug)]
22-
pub struct EpistasisTable {
23-
table: HashMap<EpistasisTableKey, HashMap<EpistasisTableKey, f64>>,
22+
pub struct EpistasisTable<S: Symbol> {
23+
table: HashMap<EpistasisTableKey<S>, HashMap<EpistasisTableKey<S>, f64>>,
2424
}
2525

26-
impl EpistasisTable {
26+
impl<S: Symbol> EpistasisTable<S> {
2727
pub fn from_model(model: &FitnessModel) -> Result<Self, VirolutionError> {
2828
let table = match &model.distribution {
2929
FitnessDistribution::Epistatic(epi_params) => {
@@ -40,16 +40,17 @@ impl EpistasisTable {
4040
}
4141

4242
pub fn from_vec(entries: Vec<EpiEntry>) -> Self {
43-
let mut table: HashMap<EpistasisTableKey, HashMap<EpistasisTableKey, f64>> = HashMap::new();
43+
let mut table: HashMap<EpistasisTableKey<S>, HashMap<EpistasisTableKey<S>, f64>> =
44+
HashMap::new();
4445
for entry in entries.iter() {
4546
table
46-
.entry((entry.pos1 as usize, Some(entry.base1)))
47+
.entry((entry.pos1 as usize, S::decode(&entry.base1)))
4748
.or_default()
48-
.insert((entry.pos2 as usize, Some(entry.base2)), entry.value);
49+
.insert((entry.pos2 as usize, S::decode(&entry.base2)), entry.value);
4950
table
50-
.entry((entry.pos2 as usize, Some(entry.base2)))
51+
.entry((entry.pos2 as usize, S::decode(&entry.base2)))
5152
.or_default()
52-
.insert((entry.pos1 as usize, Some(entry.base1)), entry.value);
53+
.insert((entry.pos1 as usize, S::decode(&entry.base1)), entry.value);
5354
}
5455
Self { table }
5556
}
@@ -68,9 +69,9 @@ impl EpistasisTable {
6869
} else {
6970
Some(EpiEntry {
7071
pos1: *pos1 as u64,
71-
base1: base1.unwrap_or_default(),
72+
base1: base1.encode(),
7273
pos2: *pos2 as u64,
73-
base2: base2.unwrap_or_default(),
74+
base2: base2.encode(),
7475
value: *value,
7576
})
7677
}
@@ -99,7 +100,7 @@ impl EpistasisTable {
99100
}
100101

101102
/// Compute factor to update the fitness of a haplotype change based on the epistasis table
102-
pub fn update_fitness(&self, haplotype: &HaplotypeRef) -> f64 {
103+
pub fn update_fitness(&self, haplotype: &HaplotypeRef<S>) -> f64 {
103104
// Callers should ensure that the haplotype is mutant
104105
if !haplotype.is_mutant() {
105106
panic!("Cannot update fitness of haplotype that is not mutant");
@@ -150,14 +151,14 @@ impl EpistasisTable {
150151
if pos != *position {
151152
if let Some(interaction) = interactions_add {
152153
// If there is an interaction, multiply the fitness
153-
if let Some(v) = interaction.get(&(*pos, *current)) {
154+
if let Some(v) = interaction.get(&(*pos, S::decode(current))) {
154155
fitness *= v;
155156
}
156157
}
157158

158159
if let Some(interaction) = interactions_remove {
159160
// If there is an interaction, divide the fitness
160-
if let Some(v) = interaction.get(&(*pos, *current)) {
161+
if let Some(v) = interaction.get(&(*pos, S::decode(current))) {
161162
fitness /= v;
162163
}
163164
}
@@ -169,16 +170,16 @@ impl EpistasisTable {
169170
}
170171

171172
/// Compute the fitness contribution within the epistasis table for a given haplotype
172-
pub fn compute_fitness(&self, haplotype: &HaplotypeRef) -> f64 {
173+
pub fn compute_fitness(&self, haplotype: &HaplotypeRef<S>) -> f64 {
173174
let mut fitness = 1.;
174175
let mutations = haplotype.get_mutations();
175176
mutations.iter().for_each(|(position, current)| {
176-
if let Some(interaction) = self.table.get(&(*position, *current)) {
177+
if let Some(interaction) = self.table.get(&(*position, S::decode(current))) {
177178
interaction
178179
.iter()
179180
.filter(|((pos, base), _)| {
180181
// Enforce that the interaction is only applied once
181-
pos > position && mutations.get(pos) != Some(base)
182+
pos > position && mutations.get(pos) != Some(&(*base.index() as u8))
182183
})
183184
.for_each(|(_, value)| {
184185
fitness *= value;

0 commit comments

Comments
 (0)