Merge branch 'add_prepmirareports' into var_of_int_format_tweak

Author: mandysulli

Cargo.toml

@@ -5,19 +5,21 @@ edition = "2024"
 description = "A set of rusty tools for use in MIRA"
 
 [dependencies]
+arrow = { version = "55.2.0", default-features = false}
+parquet = { version = "55.2.0", default-features = false, features = ["arrow"] }
 clap = { version = "4", features = ["derive"] }
 csv = "1.3.1"
 either = "1"
-serde = { version = "1.0.219", features = ["derive"] }
-serde_yaml = "0.9"
 glob = "0.3.2"
 ordered-float = "5.0.0"
 #plotly = "0.12.1"
 plotly = { git = "https://github.com/plotly/plotly.rs.git", branch = "main" }
-
+serde = { version = "1.0.219", features = ["derive"] }
+serde_json = "1.0"
+serde_yaml_ng = "0.10.0"
 zoe = { version = "0.0.19", default-features = false, features = [
     "multiversion",
 ] }
 
 [profile.release]
-strip = true
+strip = true

src/main.rs

@@ -28,6 +28,8 @@ enum Commands {
     NTDiffs(NTDiffsArgs),
     /// Plotter
     Plotter(PlotterArgs),
+    /// Prepare MIRA report
+    PrepareMiraReports(ReportsArgs),
 }
 
 fn main() {
@@ -40,22 +42,25 @@ fn main() {
         }
         Commands::PositionsOfInterest(cmd_args) => positions_of_interest_process(cmd_args)
             .expect(&format!("{module}::PositionsOfInterest")),
+
         Commands::FindChemistry(cmd_args) => {
             find_chemistry_process(cmd_args).unwrap_or_die(&format!("{module}::FindChemistry"))
         }
         Commands::Hamming(cmd_args) => {
-            all_sample_hd_process(cmd_args).unwrap_or_die(&format!("{module}::Hamming"))
+            all_sample_hd_process(cmd_args).unwrap_or_die(&format!("{module}::Hamming"));
         }
         Commands::NTDiffs(cmd_args) => all_sample_nt_diffs_process(cmd_args),
         Commands::Plotter(cmd_args) => {
-            plotter_process(cmd_args).expect(&format!("{module}::Plotter"))
+            plotter_process(cmd_args).unwrap_or_else(|_| panic!("{module}::Plotter"))
         }
+        Commands::PrepareMiraReports(cmd_args) => prepare_mira_reports_process(cmd_args)
+            .unwrap_or_else(|_| panic!("{module}::PrepareMiraReports")),
         _ => {
             eprintln!("mira-oxide: unrecognized command {:?}", args.command);
             std::process::exit(1)
         }
     }
 }
 
-mod processes;
-pub use crate::processes::*;
+pub mod processes;
+pub mod utils;

src/processes/all_sample_nt_diffs.rs

@@ -5,10 +5,7 @@ use std::{
     io::{BufReader, BufWriter, Write, stdin, stdout},
     path::PathBuf,
 };
-use zoe::{
-    data::fasta::FastaNT,
-    prelude::*,
-};
+use zoe::{data::fasta::FastaNT, prelude::*};
 
 #[derive(Debug, Parser)]
 #[command(
@@ -70,39 +67,36 @@ pub fn all_sample_nt_diffs_process(args: NTDiffsArgs) {
             record.map(|r| {
                 let FastaNT { name, sequence } = r.recode_to_dna();
                 ValidSeq {
-                    name, 
+                    name,
                     sequence,
                 }
               }))
         .collect::<Result<Vec<_>, _>>()
         .unwrap_or_die("Could not process other data.");
-    
-        writeln!(
-            &mut writer,
-            "sequence_1{}sequence_2{}nt_sequence_1{}position{}nt_sequence_2",
-            delim, delim, delim, delim
-        ).unwrap();
 
+    writeln!(
+        &mut writer,
+        "sequence_1{delim}sequence_2{delim}nt_sequence_1{delim}position{delim}nt_sequence_2"
+    )
+    .unwrap();
+
+    all_sequences.iter().for_each(|f| {
+        let name_1 = &f.name;
+        let seq1 = &f.sequence;
         all_sequences.iter().for_each(|f| {
-            let name_1 = &f.name;
-            let seq1 = &f.sequence;
-            all_sequences.iter().for_each(|f| {
-                let name_2 = &f.name;
-                let seq2 = &f.sequence;
-                for (i, (nt1, nt2)) in seq1.iter().zip(seq2.iter()).enumerate() {
-                    if nt1 != nt2 {
-                        let nucleotide1 = char::from(*nt1);
-                        let nucleotide2 = char::from(*nt2);
-                        writeln!(
-                            &mut writer,
-                            "{}{}{}{}{}{}{}{}{}",
-                            name_1, delim, name_2, delim, nucleotide1, delim, i, delim, nucleotide2
-                        )
-                        .unwrap();
-                    }
+            let name_2 = &f.name;
+            let seq2 = &f.sequence;
+            for (i, (nt1, nt2)) in seq1.iter().zip(seq2.iter()).enumerate() {
+                if nt1 != nt2 {
+                    let nucleotide1 = char::from(*nt1);
+                    let nucleotide2 = char::from(*nt2);
+                    writeln!(
+                        &mut writer,
+                        "{name_1}{delim}{name_2}{delim}{nucleotide1}{delim}{i}{delim}{nucleotide2}"
+                    )
+                    .unwrap();
                 }
-            });
+            }
         });
-        
-
+    });
 }

src/processes/find_chemistry.rs

@@ -228,7 +228,7 @@ fn get_config_path(args: &FindChemArgs, seq_len: Option<usize>) -> String {
         .wd_path
         .to_str()
         .expect("Failed to convert work directory path to string");
-    format!("{}{}", wd_path, path_extension)
+    format!("{wd_path}{path_extension}")
 }
 
 #[derive(Debug)]
@@ -274,7 +274,7 @@ impl fmt::Display for ChemistryOutput {
 /// sequences, returns None
 fn get_average_line_length(fastq: &PathBuf) -> Result<Option<usize>, std::io::Error> {
     let sample_size = 5;
-    let file = File::open(&fastq)?;
+    let file = File::open(fastq)?;
     let buf_reader = BufReader::new(file);
     let fastq_reader = FastQReader::new(buf_reader);
 
@@ -313,7 +313,7 @@ pub fn find_chemistry_process(args: FindChemArgs) -> Result<(), std::io::Error>
     //let args = CheckChemArgs::parse();
     // handle input validation to ensure valid combinations of
     if let Err(e) = args.validate() {
-        eprintln!("Error: {}", e);
+        eprintln!("Error: {e}");
         std::process::exit(1);
     }
     // parse the arguments into output format

src/processes/mod.rs

@@ -3,4 +3,5 @@ pub mod all_sample_nt_diffs;
 pub mod find_chemistry;
 pub mod plotter;
 pub mod positions_of_interest;
+pub mod prepare_mira_reports;
 pub mod variants_of_interest;

src/processes/plotter.rs

@@ -162,7 +162,7 @@ fn generate_plot_coverage(input_directory: &PathBuf) -> Result<Plot, Box<dyn Err
 
     // Set the figure title
     let layout = Layout::new()
-        .title(&format!(
+        .title(format!(
             "Coverage | {}",
             input_directory
                 .file_name()
@@ -201,61 +201,59 @@ fn generate_plot_coverage_seg(input_directory: &PathBuf) -> Result<Plot, Box<dyn
     let mut file_paths = Vec::new();
 
     // First, count files and collect paths
-    for entry in glob(&format!(
+    for path in (glob(&format!(
         "{}/tables/*coverage.txt",
         input_directory.display()
-    ))? {
-        if let Ok(path) = entry {
-            //file_count += 1;
-            file_paths.push(path);
-        }
+    ))?)
+    .flatten()
+    {
+        //file_count += 1;
+        file_paths.push(path);
     }
 
     // Calculate grid dimensions for subplots
     let rows = 4; //((file_count + 2) as f64).sqrt().ceil() as usize;
     let cols = 2; //(file_count + rows - 1) / rows; // Ceiling division
 
     // Load variant data into a HashMap keyed by segment name
+    // TODO: consider a struct with named fields
     let mut variants_data: HashMap<String, Vec<(u32, String, String, u32, u32, f32)>> =
         HashMap::new();
 
     // Look for variant files with matching prefixes in the directory
-    for entry in glob(&format!(
+    for variant_path in (glob(&format!(
         "{}/tables/*variants.txt",
         input_directory.display()
-    ))? {
-        if let Ok(variant_path) = entry {
-            let file = File::open(&variant_path)?;
-
-            // Create a TSV reader
-            let mut rdr = ReaderBuilder::new()
-                .delimiter(b'\t')
-                .has_headers(true)
-                .from_reader(file);
-
-            for result in rdr.records() {
-                let record = result?;
-                if record.len() >= 8 {
-                    let segment_name = record[0].to_string();
-                    let position: u32 = record[1].parse()?;
-                    let consensus_allele: String = record[3].to_string();
-                    let minority_allele: String = record[4].to_string();
-                    let consensus_count: u32 = record[5].parse()?;
-                    let minority_count: u32 = record[6].parse()?;
-                    let minority_frequency: f32 = record[8].parse()?;
-
-                    variants_data
-                        .entry(segment_name)
-                        .or_insert_with(Vec::new)
-                        .push((
-                            position,
-                            consensus_allele,
-                            minority_allele,
-                            consensus_count,
-                            minority_count,
-                            minority_frequency,
-                        ));
-                }
+    ))?)
+    .flatten()
+    {
+        let file = File::open(&variant_path)?;
+
+        // Create a TSV reader
+        let mut rdr = ReaderBuilder::new()
+            .delimiter(b'\t')
+            .has_headers(true)
+            .from_reader(file);
+
+        for result in rdr.records() {
+            let record = result?;
+            if record.len() >= 8 {
+                let segment_name = record[0].to_string();
+                let position: u32 = record[1].parse()?;
+                let consensus_allele: String = record[3].to_string();
+                let minority_allele: String = record[4].to_string();
+                let consensus_count: u32 = record[5].parse()?;
+                let minority_count: u32 = record[6].parse()?;
+                let minority_frequency: f32 = record[8].parse()?;
+
+                variants_data.entry(segment_name).or_default().push((
+                    position,
+                    consensus_allele,
+                    minority_allele,
+                    consensus_count,
+                    minority_count,
+                    minority_frequency,
+                ));
             }
         }
     }
@@ -356,7 +354,7 @@ fn generate_plot_coverage_seg(input_directory: &PathBuf) -> Result<Plot, Box<dyn
             // Create trace for minority values with consistent color (but with transparency)
             let minority_trace = Scatter::new(variant_positions, minority_values)
                 .mode(Mode::Markers)
-                .name(&format!("{}", segment_name))
+                .name(&segment_name)
                 .marker(
                     plotly::common::Marker::new()
                         .color(segment_color)
@@ -383,7 +381,7 @@ fn generate_plot_coverage_seg(input_directory: &PathBuf) -> Result<Plot, Box<dyn
                 .columns(cols)
                 .pattern(GridPattern::Independent),
         )
-        .title(&format!(
+        .title(format!(
             "Segment Coverage | {}",
             input_directory
                 .file_name()
@@ -505,17 +503,15 @@ fn generate_sankey_plot(input_directory: &PathBuf) -> Result<Plot, Box<dyn Error
 
     // Process data and build node map first
     let mut records = Vec::new();
-    for line in lines {
-        if let Ok(line) = line {
-            let parts: Vec<&str> = line.split('\t').collect();
-            if parts.len() >= 3 {
-                let record = parts[0];
-                let reads: u32 = parts[1].parse().unwrap_or(0);
-
-                // Skip "NA" values and 0 reads
-                if parts[1] != "NA" && reads > 0 {
-                    records.push((record.to_string(), reads));
-                }
+    for line in lines.map_while(Result::ok) {
+        let parts: Vec<&str> = line.split('\t').collect();
+        if parts.len() >= 3 {
+            let record = parts[0];
+            let reads: u32 = parts[1].parse().unwrap_or(0);
+
+            // Skip "NA" values and 0 reads
+            if parts[1] != "NA" && reads > 0 {
+                records.push((record.to_string(), reads));
             }
         }
     }
@@ -733,7 +729,7 @@ fn generate_sankey_plot(input_directory: &PathBuf) -> Result<Plot, Box<dyn Error
 
     // Set layout
     let layout = Layout::new()
-        .title(&format!(
+        .title(format!(
             "Read Assignment | {}",
             input_directory
                 .file_name()