dev/fixes (#32)

* fix(kraken): include unclassified reads in the tree building step

* feat: allow missing taxids to be supplied

* feat: early bail if no reads found

* feat: add missing taxons to the final report

fixes clippy warnings for complex types by introducing structs for processed kraken outputs and processed kraken trees as well

* docs(changelog): update changelog for 3.0.0

* refactor: reorder report format

* docs(readme): update for changes

Author: Sam-Sims

CHANGELOG.md

@@ -5,6 +5,18 @@ All notable changes to this project will be documented in this file.
 The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.1.0/),
 and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
 
+## [3.0.0] - 2025-09-26
+
+### Added
+- Able to specifiy taxon ids that are not present, without kractor stopping. These are instead logged to stderr with a warning. This may be useful when running kractor in a wrapper script for several fastq files and just want to extract a set of taxonids from them all - without caring if they are present or not.
+- Include a new field `missing_taxon_ids` in the summary output.
+
+### Changed
+- Under the hood refactoring, introducing structs for the processed kraken outputs and processed kraken trees to simplify the returned data.
+
+### Fixed
+- Unclassified reads being skipped in the tree building stage, meaning they were unable to be extracted 
+
 ## [2.0.0] - 2025-08-12
 
 ### Added
@@ -102,4 +114,4 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 ## [0.1.0]
 
 ### Added
-- Initial release
+- Initial release

README.md

@@ -150,6 +150,9 @@ Use `--summary` to get summary statistics (output to stdout on completion)
 ```json
 {
   "total_taxon_count": 2,
+  "missing_taxon_ids": [
+    999999999
+  ]
   "reads_extracted_per_taxon": {
     "0": 745591,
     "1": 1646
@@ -199,7 +202,9 @@ One or more taxonomic IDs to extract.
 
 For example: `-t 1 2 10`
 
-Each taxid is affected by `--exclude`, `--parents`, and `--children` if those options are used.
+Each taxonomic id is affected by `--exclude`, `--parents`, and `--children` if those options are used.
+
+Taxonomic ids do not need to be present in a given report. This may be useful when running kractor in a wrapper script for several fastq files and just want to extract a set of taxon ids from them all - without caring if they are present or not.
 
 ### Optional:
 

src/extract.rs

@@ -1,12 +1,20 @@
 use crate::parsers::fastx::{parse_fastq, write_output_fasta, write_output_fastq};
-use crate::parsers::kraken::{build_tree_from_kraken_report, extract_children, extract_parents};
+use crate::parsers::kraken::{
+    build_tree_from_kraken_report, extract_children, extract_parents, ProcessedKrakenTree,
+};
 use color_eyre::{eyre::bail, eyre::eyre, eyre::WrapErr, Result};
 use crossbeam::{channel, thread};
 use fxhash::FxHashSet;
-use log::debug;
+use log::{debug, info, warn};
 use noodles::fastq;
 use std::path::PathBuf;
 
+#[derive(Debug, Clone)]
+pub struct CollectedTaxonIds {
+    pub found: Vec<i32>,
+    pub missing: Vec<i32>,
+}
+
 pub fn process_single_end(
     reads_to_save: &FxHashSet<Vec<u8>>,
     input: &[PathBuf],
@@ -139,39 +147,56 @@ pub fn collect_taxons_to_save(
     children: bool,
     parents: bool,
     taxids: &[i32],
-) -> Result<Vec<i32>> {
+) -> Result<CollectedTaxonIds> {
     let mut taxon_ids_to_save = Vec::new();
+    let mut missing_taxon_ids = Vec::new();
 
     // I dont think we will reach this code ever since clap should catch this - but in case it doesnt
     if (parents || children) && report.is_none() {
         return Err(eyre!("Report required when parents or children is enabled"));
     }
 
     if let Some(report_path) = report {
-        let (nodes, taxon_map) = build_tree_from_kraken_report(taxids, report_path)
+        let ProcessedKrakenTree {
+            nodes,
+            taxon_map,
+            missing_taxon_ids: missing_ids,
+        } = build_tree_from_kraken_report(taxids, report_path)
             .wrap_err("Failed to build tree from Kraken report")?;
 
+        if !missing_ids.is_empty() {
+            warn!(
+                "The following taxon IDs were not found in the kraken report and will be ignored: {:?}",
+                missing_ids
+            );
+        }
+
+        missing_taxon_ids = missing_ids;
+
+        // remove missing taxon ids from the input list
+        let taxids: Vec<i32> = taxids
+            .iter()
+            .filter(|id| !missing_taxon_ids.contains(id))
+            .cloned()
+            .collect();
+
+        if taxon_map.is_empty() {
+            bail!("No valid taxon IDs found in the kraken report");
+        }
+
         if children {
             debug!("Extracting children");
             let mut children = Vec::new();
             for taxid in taxids {
-                if let Some(&node_index) = taxon_map.get(taxid) {
-                    extract_children(&nodes, node_index, &mut children).wrap_err_with(|| {
-                        format!("Failed to extract children for taxon ID {taxid}")
-                    })?;
-                } else {
-                    return Err(eyre!("Taxon ID {} not found in taxonomy map", taxid));
+                if let Some(&node_index) = taxon_map.get(&taxid) {
+                    extract_children(&nodes, node_index, &mut children)?;
                 }
             }
             taxon_ids_to_save.extend(children);
         } else if parents {
             debug!("Extracting parents");
             for taxid in taxids {
-                taxon_ids_to_save.extend(
-                    extract_parents(&taxon_map, &nodes, *taxid).wrap_err_with(|| {
-                        format!("Failed to extract parents for taxon ID {taxid}")
-                    })?,
-                );
+                taxon_ids_to_save.extend(extract_parents(&taxon_map, &nodes, taxid)?);
             }
         } else {
             taxon_ids_to_save.extend(taxids);
@@ -184,10 +209,18 @@ pub fn collect_taxons_to_save(
     taxon_ids_to_save.sort_unstable();
     taxon_ids_to_save.dedup();
 
+    missing_taxon_ids.sort_unstable();
+    missing_taxon_ids.dedup();
+
     if taxon_ids_to_save.is_empty() {
         bail!("No taxon IDs were identified for extraction");
     }
-    Ok(taxon_ids_to_save)
+
+    info!("Identified {} taxon IDs to save", taxon_ids_to_save.len());
+    Ok(CollectedTaxonIds {
+        found: taxon_ids_to_save,
+        missing: missing_taxon_ids,
+    })
 }
 
 #[cfg(test)]
@@ -383,48 +416,58 @@ mod tests {
     #[test]
     fn test_no_report() {
         let taxids = vec![123, 456, 789];
-        let saved_taxons = collect_taxons_to_save(&None, false, false, &taxids).unwrap();
+        let collected = collect_taxons_to_save(&None, false, false, &taxids).unwrap();
+
+        assert_eq!(collected.found, taxids);
+        assert!(collected.missing.is_empty());
+    }
+
+    #[test]
+    fn test_extract_unclassified() {
+        let dir = tempdir().unwrap();
+        let report_path = create_test_kraken_report(&dir);
+        let taxids = vec![0, 2];
+        let collected = collect_taxons_to_save(&Some(report_path), false, false, &taxids).unwrap();
 
-        assert_eq!(saved_taxons, taxids);
+        assert_eq!(collected.found, vec![0, 2]);
+        assert!(collected.missing.is_empty());
     }
 
     #[test]
     fn test_no_parents_no_children() {
         let dir = tempdir().unwrap();
         let report_path = create_test_kraken_report(&dir);
         let taxids = vec![1385, 1386, 91061];
-        let saved_taxons =
-            collect_taxons_to_save(&Some(report_path), false, false, &taxids).unwrap();
+        let collected = collect_taxons_to_save(&Some(report_path), false, false, &taxids).unwrap();
 
-        assert_eq!(saved_taxons, taxids);
+        assert_eq!(collected.found, taxids);
+        assert!(collected.missing.is_empty());
     }
 
     #[test]
     fn test_children() {
         let dir = tempdir().unwrap();
         let report_path = create_test_kraken_report(&dir);
         let taxids = vec![1239];
-        let saved_taxons =
-            collect_taxons_to_save(&Some(report_path), true, false, &taxids).unwrap();
+        let collected = collect_taxons_to_save(&Some(report_path), true, false, &taxids).unwrap();
 
-        assert!(saved_taxons.contains(&1239));
-        assert!(saved_taxons.contains(&91062));
-        assert!(saved_taxons.contains(&91061));
+        assert!(collected.found.contains(&1239));
+        assert!(collected.found.contains(&91062));
+        assert!(collected.found.contains(&91061));
     }
 
     #[test]
     fn test_parents() {
         let dir = tempdir().unwrap();
         let report_path = create_test_kraken_report(&dir);
         let taxids = vec![91061];
-        let saved_taxons =
-            collect_taxons_to_save(&Some(report_path), false, true, &taxids).unwrap();
-
-        assert!(saved_taxons.contains(&91061));
-        assert!(saved_taxons.contains(&1239));
-        assert!(saved_taxons.contains(&1783272));
-        assert!(saved_taxons.contains(&131567));
-        assert!(saved_taxons.contains(&2));
+        let collected = collect_taxons_to_save(&Some(report_path), false, true, &taxids).unwrap();
+
+        assert!(collected.found.contains(&91061));
+        assert!(collected.found.contains(&1239));
+        assert!(collected.found.contains(&1783272));
+        assert!(collected.found.contains(&131567));
+        assert!(collected.found.contains(&2));
     }
 
     #[test]
@@ -437,12 +480,24 @@ mod tests {
         assert!(result.is_err());
     }
 
+    #[test]
+    fn test_missing_taxons_recorded() {
+        let dir = tempdir().unwrap();
+        let report_path = create_test_kraken_report(&dir);
+        let taxids = vec![1239, 999];
+        let collected = collect_taxons_to_save(&Some(report_path), false, false, &taxids).unwrap();
+
+        assert!(collected.found.contains(&1239));
+        assert_eq!(collected.missing, vec![999]);
+    }
+
     #[test]
     fn test_dedup_and_sort() {
         let taxids = vec![456, 123, 456, 789, 123];
-        let saved_taxons = collect_taxons_to_save(&None, false, false, &taxids).unwrap();
+        let collected = collect_taxons_to_save(&None, false, false, &taxids).unwrap();
 
-        assert_eq!(saved_taxons, vec![123, 456, 789]);
+        assert_eq!(collected.found, vec![123, 456, 789]);
+        assert!(collected.missing.is_empty());
     }
 
     #[test]

src/kractor.rs

@@ -1,9 +1,10 @@
 use crate::extract::{process_paired_end, process_single_end};
+use crate::parsers::kraken::ProcessedKrakenOutput;
 use crate::{extract, parsers, Cli};
-use color_eyre::eyre::ensure;
+use color_eyre::eyre::{bail, ensure};
 use color_eyre::Result;
 use fxhash::{FxHashMap, FxHashSet};
-use log::{debug, info};
+use log::info;
 use serde::{Deserialize, Serialize};
 
 #[derive(Serialize, Deserialize)]
@@ -16,11 +17,13 @@ struct Summary {
     input_format: String,
     output_format: String,
     kractor_version: String,
+    missing_taxon_ids: Vec<i32>,
 }
 
 pub struct Kractor {
     args: Cli,
     taxon_ids: Vec<i32>,
+    missing_taxon_ids: Vec<i32>,
     reads_to_save: FxHashSet<Vec<u8>>,
     reads_per_taxon: FxHashMap<i32, usize>,
     summary: Option<Summary>,
@@ -31,6 +34,7 @@ impl Kractor {
         Self {
             args,
             taxon_ids: Vec::new(),
+            missing_taxon_ids: Vec::new(),
             reads_to_save: FxHashSet::default(),
             reads_per_taxon: FxHashMap::default(),
             summary: None,
@@ -49,24 +53,34 @@ impl Kractor {
     }
 
     fn collect_taxons(&mut self) -> Result<()> {
-        self.taxon_ids = extract::collect_taxons_to_save(
+        let collected = extract::collect_taxons_to_save(
             &self.args.report,
             self.args.children,
             self.args.parents,
             &self.args.taxid,
         )?;
-        debug!("Taxon IDs identified: {:?}", self.taxon_ids);
+        self.taxon_ids = collected.found;
+        self.missing_taxon_ids = collected.missing;
         Ok(())
     }
 
     fn process_kraken_output(&mut self) -> Result<()> {
-        (self.reads_to_save, self.reads_per_taxon) = parsers::kraken::process_kraken_output(
+        let ProcessedKrakenOutput {
+            reads_to_save,
+            reads_per_taxon,
+        } = parsers::kraken::process_kraken_output(
             &self.args.kraken,
             self.args.exclude,
             &self.taxon_ids,
         )?;
+        self.reads_to_save = reads_to_save;
+        self.reads_per_taxon = reads_per_taxon;
 
-        debug!("Identified {} reads to save", self.reads_to_save.len());
+        if self.reads_to_save.is_empty() {
+            bail!("No reads found for the specified taxon ID(s). Nothing to extract.");
+        }
+
+        info!("Identified {} reads to save", self.reads_to_save.len());
         Ok(())
     }
 
@@ -102,6 +116,7 @@ impl Kractor {
                     "fastq".to_string()
                 },
                 kractor_version: env!("CARGO_PKG_VERSION").to_string(),
+                missing_taxon_ids: self.missing_taxon_ids.clone(),
             });
         } else {
             let (reads_parsed1, reads_output1) = process_single_end(
@@ -119,6 +134,7 @@ impl Kractor {
             self.summary = Some(Summary {
                 total_taxon_count: self.taxon_ids.len(),
                 reads_extracted_per_taxon: self.reads_per_taxon.clone(),
+                missing_taxon_ids: self.missing_taxon_ids.clone(),
                 total_reads_in: reads_in,
                 total_reads_out: reads_out,
                 proportion_extracted: reads_out as f64 / reads_in as f64,
@@ -152,7 +168,6 @@ impl Kractor {
         );
         self.validate_outputs()?;
         self.collect_taxons()?;
-        info!("{} taxons identified to save", self.taxon_ids.len());
         info!("Processing Kraken2 output file");
         self.process_kraken_output()?;
         info!("Processing reads");