with --debug, store time information regarding time taken to genotype…

… repeat

this requires repeat instances to be mutable
this is done to be able to identify repeats that are slow to genotype
may revert later

Cargo.toml

@@ -1,6 +1,6 @@
 [package]
 name = "STRdust"
-version = "0.8.6"
+version = "0.8.7"
 edition = "2021"
 
 
@@ -28,6 +28,7 @@ petgraph = "0.6.4"
 hts-sys = "2.1.1"
 reqwest = { version = "0.12.5", features = ["blocking", "json"] }
 indicatif = { version = "0.17.1", features = ["rayon"] }
+chrono = "0.4.38"
 
 [dev-dependencies]
 ctor = "*"

src/call.rs

@@ -21,8 +21,8 @@ pub fn genotype_repeats(args: Cli) {
         // The indexedreader is created once and passed on to the function
         let num_intervals = repeats.num_intervals;
         let mut bam = parse_bam::create_bam_reader(&args.bam, &args.fasta);
-        for repeat in repeats.progress_count(num_intervals as u64) {
-            if let Ok(output) = genotype::genotype_repeat_singlethreaded(&repeat, &args, &mut bam) {
+        for mut repeat in repeats.progress_count(num_intervals as u64) {
+            if let Ok(output) = genotype::genotype_repeat_singlethreaded(&mut repeat, &args, &mut bam) {
                 writeln!(handle, "{output}").expect("Failed writing the result.");
             }
         }
@@ -39,8 +39,8 @@ pub fn genotype_repeats(args: Cli) {
         repeats
             .par_bridge()
             .progress_count(num_intervals as u64)
-            .for_each(|repeat| {
-                if let Ok(output) = genotype::genotype_repeat_multithreaded(&repeat, &args) {
+            .for_each(|mut repeat| {
+                if let Ok(output) = genotype::genotype_repeat_multithreaded(&mut repeat, &args) {
                     let mut geno = genotypes.lock().expect("Unable to lock genotypes mutex");
                     geno.push(output);
                 } else {
@@ -66,8 +66,8 @@ pub fn genotype_repeats(args: Cli) {
         repeats
             .par_bridge()
             .progress_count(num_intervals as u64)
-            .for_each(|repeat| {
-                if let Ok(output) = genotype::genotype_repeat_multithreaded(&repeat, &args) {
+            .for_each(|mut repeat| {
+                if let Ok(output) = genotype::genotype_repeat_multithreaded(&mut repeat, &args) {
                     println!("{output}");
                 }
             });

src/consensus.rs

@@ -197,6 +197,7 @@ mod tests {
                 chrom: "chr1".to_string(),
                 start: 1,
                 end: 100,
+                created: None
             },
         );
         println!("Consensus: {}", cons.seq.unwrap());

src/genotype.rs

@@ -7,19 +7,25 @@ use rust_htslib::bam;
 // when running multithreaded, the indexedreader has to be created every time again
 // this is probably expensive
 pub fn genotype_repeat_multithreaded(
-    repeat: &crate::repeats::RepeatInterval,
+    repeat: &mut crate::repeats::RepeatInterval,
     args: &Cli,
 ) -> Result<crate::vcf::VCFRecord, String> {
+    if args.debug {
+        repeat.set_time_stamp()
+    }
     let mut bam = parse_bam::create_bam_reader(&args.bam, &args.fasta);
     genotype_repeat(repeat, args, &mut bam)
 }
 
 // when running singlethreaded, the indexedreader is created once and simply passed on
 pub fn genotype_repeat_singlethreaded(
-    repeat: &crate::repeats::RepeatInterval,
+    repeat: &mut crate::repeats::RepeatInterval,
     args: &Cli,
     bam: &mut bam::IndexedReader,
 ) -> Result<crate::vcf::VCFRecord, String> {
+    if args.debug {
+        repeat.set_time_stamp()
+    }
     genotype_repeat(repeat, args, bam)
 }
 
@@ -42,18 +48,19 @@ fn genotype_repeat(
                 repeat,
                 "N",
                 ".".to_string(),
+                args
             ))
         }
     };
 
     let repeat_compressed_reference = repeat.make_repeat_compressed_sequence(&args.fasta, flanking);
-    if args.debug {
-        // write the repeat compressed reference to a file
-        use std::fs;
-        let header = format!(">{chrom}:{start}-{end}\n", chrom = repeat.chrom, start = repeat.start, end = repeat.end);
-        let fas = String::from_utf8(repeat_compressed_reference.clone()).expect("Unable to convert repeat compressed reference to string");
-        fs::write("repeat_compressed.fa", header + &fas).expect("Unable to write repeat compressed reference to file");
-    }
+    // if args.debug {
+    //     // write the repeat compressed reference to a file
+    //     use std::fs;
+    //     let header = format!(">{chrom}:{start}-{end}\n", chrom = repeat.chrom, start = repeat.start, end = repeat.end);
+    //     let fas = String::from_utf8(repeat_compressed_reference.clone()).expect("Unable to convert repeat compressed reference to string");
+    //     fs::write("repeat_compressed.fa", header + &fas).expect("Unable to write repeat compressed reference to file");
+    // }
 
     // alignments can be extracted in an unphased manner, if the chromosome is --haploid or the --unphased is set
     // this means that --haploid overrides the phases which could be present in the bam file
@@ -69,6 +76,7 @@ fn genotype_repeat(
                 repeat,
                 &repeat_ref_seq,
                 0.to_string(),
+                args
             ));
         }
     };
@@ -119,6 +127,7 @@ fn genotype_repeat(
                 repeat,
                 &repeat_ref_seq,
                 insertions.len().to_string(),
+                args
             ));
         }
         // there is only one haplotype, haploid, so this gets duplicated for reporting in the VCF module
@@ -144,6 +153,7 @@ fn genotype_repeat(
                 repeat,
                 &repeat_ref_seq,
                 insertions.len().to_string(),
+                args
             ));
         }
         // handle the special case where there is only 1 read for the repeat (but minimal support is 1 so it passes)
@@ -156,6 +166,7 @@ fn genotype_repeat(
                 all_insertions,
                 reads.ps,
                 flags,
+                args
             ));
         }
 
@@ -231,6 +242,7 @@ fn genotype_repeat(
         repeat,
         reads.ps,
         flags,
+        args
     ))
 }
 
@@ -344,6 +356,7 @@ mod tests {
             chrom: String::from("chr7"),
             start: 154654404,
             end: 154654432,
+            created: None
         };
         let flanking = 2000;
         let minlen = 5;
@@ -383,6 +396,7 @@ mod tests {
             chrom: String::from("chr7"),
             start: 154654404,
             end: 154654432,
+            created: None
         };
         let args = Cli {
             bam: String::from("test_data/small-test-phased.bam"),
@@ -412,6 +426,7 @@ mod tests {
             chrom: String::from("chr7"),
             start: 154654404,
             end: 154654432,
+            created: None
         };
         let args = Cli {
             bam: String::from("test_data/small-test-phased.bam"),
@@ -458,6 +473,7 @@ mod tests {
             chrom: String::from("chr7"),
             start: 154654404,
             end: 154654432,
+            created: None
         };
         let mut bam = parse_bam::create_bam_reader(&args.bam, &args.fasta);
         let genotype = genotype_repeat(&repeat, &args, &mut bam);
@@ -487,6 +503,7 @@ mod tests {
             chrom: String::from("chr7"),
             start: 154654404,
             end: 154654432,
+            created: None
         };
         let mut bam = parse_bam::create_bam_reader(&args.bam, &args.fasta);
         let genotype = genotype_repeat(&repeat, &args, &mut bam);
@@ -517,6 +534,7 @@ mod tests {
             chrom: String::from("chr7"),
             start: 154654404,
             end: 154654432,
+            created: None
         };
         let mut bam = parse_bam::create_bam_reader(&args.bam, &args.fasta);
         let genotype = genotype_repeat(&repeat, &args, &mut bam);

src/parse_bam.rs

@@ -124,6 +124,7 @@ fn test_get_overlapping_reads() {
         chrom: String::from("chr7"),
         start: 154654404,
         end: 154654432,
+        created: None,
     };
     let unphased = false;
     let mut bam = create_bam_reader(&bam, &fasta);
@@ -138,6 +139,7 @@ fn test_get_overlapping_reads_url() {
         chrom: String::from("chr20"),
         start: 154654404,
         end: 154654432,
+        created: None
     };
     let unphased = false;
     let mut bam = create_bam_reader(&bam, &fasta);

src/phase_insertions.rs

@@ -332,6 +332,7 @@ mod tests {
                 chrom: "chr7".to_string(),
                 start: 154654404,
                 end: 154654432,
+                created: None
             },
             false,
         );
@@ -372,6 +373,7 @@ mod tests {
                 chrom: "chr7".to_string(),
                 start: 154654404,
                 end: 154654432,
+                created: None
             },
             false,
         );
@@ -411,6 +413,7 @@ mod tests {
                 chrom: "chr7".to_string(),
                 start: 154654404,
                 end: 154654432,
+                created: None
             },
             false,
         );
@@ -453,6 +456,7 @@ mod tests {
                 chrom: "chr7".to_string(),
                 start: 154654404,
                 end: 154654432,
+                created: None
             },
             false,
         );
@@ -493,6 +497,7 @@ mod tests {
                 chrom: "chr7".to_string(),
                 start: 154654404,
                 end: 154654432,
+                created: None
             },
             false,
         );
@@ -564,6 +569,7 @@ mod tests {
                 chrom: "chr7".to_string(),
                 start: 154654404,
                 end: 154654432,
+                created: None
             },
             false,
         );

src/repeats.rs

@@ -80,6 +80,7 @@ impl Clone for RepeatInterval {
             chrom: self.chrom.clone(),
             start: self.start,
             end: self.end,
+            created: self.created,
         }
     }
 }
@@ -105,6 +106,7 @@ pub struct RepeatInterval {
     pub chrom: String,
     pub start: u32,
     pub end: u32,
+    pub created: Option<chrono::DateTime<chrono::Utc>>,
 }
 
 impl fmt::Display for RepeatInterval {
@@ -144,7 +146,7 @@ impl RepeatInterval {
                     .expect("Failed parsing chromosome length from fai file")
                     > end
             {
-                return Some(Self { chrom, start, end });
+                return Some(Self { chrom, start, end, created: None });
             }
         }
         // if the chromosome is not in the fai file or the end does not fit the interval, return None
@@ -157,6 +159,7 @@ impl RepeatInterval {
             chrom: chrom.to_string(),
             start,
             end,
+            created: None
         }
     }
 
@@ -201,6 +204,10 @@ impl RepeatInterval {
         }
         Some(repeat_ref_sequence)
     }
+
+    pub fn set_time_stamp(& mut self) {
+        self.created = Some(chrono::Utc::now());
+    }
 }
 
 #[cfg(test)]