diff --git a/src/bin/coverm.rs b/src/bin/coverm.rs index 91a0070..70769c5 100644 --- a/src/bin/coverm.rs +++ b/src/bin/coverm.rs @@ -1679,16 +1679,32 @@ impl FilterParameters { pub fn add_metabat_filtering_if_required(&mut self, m: &clap::ArgMatches) { if doing_metabat(m) { + // A threshold given with --min-read-percent-identity is kept; the + // 97% default only replaces an unset one. + if m.contains_id("min-read-percent-identity") { + info!( + "Using the given single read percent identity threshold of {} for \ + MetaBAT adjusted coverage", + self.min_percent_identity_single + ); + } else { + // we use >= where metabat uses >. Gah. + self.min_percent_identity_single = 0.97001; + info!( + "Setting single read percent identity threshold at 0.97 for \ + MetaBAT adjusted coverage" + ); + } + // Improper pair and supplementary alignments are already included + // unless the user asked otherwise with --proper-pairs-only or + // --exclude-supplementary, which are respected. Secondary + // alignments are excluded by default, so include them here. + self.flag_filters.include_secondary = true; info!( - "Setting single read percent identity threshold at 0.97 for \ - MetaBAT adjusted coverage, and not filtering out supplementary, \ - secondary and improper pair alignments" + "For MetaBAT adjusted coverage, including secondary alignments; \ + supplementary alignments included: {}; improper pair alignments included: {}", + self.flag_filters.include_supplementary, self.flag_filters.include_improper_pairs ); - // we use >= where metabat uses >. Gah. - self.min_percent_identity_single = 0.97001; - self.flag_filters.include_improper_pairs = true; - self.flag_filters.include_supplementary = true; - self.flag_filters.include_secondary = true; } } diff --git a/src/cli.rs b/src/cli.rs index 5d6571d..4b15123 100644 --- a/src/cli.rs +++ b/src/cli.rs @@ -208,8 +208,15 @@ fn add_thresholding_options(manual: Manual) -> Manual { .long("--min-read-percent-identity") .help(&format!( "Exclude reads by overall percent \ - identity e.g. 95 for 95%. {}", - default_roff("0") + identity e.g. 95 for 95%. Reads at exactly the threshold \ + are kept ({}). {} With {} in {} the default is {} \ + (see the {} method); a given value is used as is.", + &monospace_roff(">="), + default_roff("0"), + &monospace_roff("--methods metabat"), + &monospace_roff("coverm contig"), + &monospace_roff("97"), + &monospace_roff("metabat") )), ) .option( @@ -965,7 +972,13 @@ pub fn contig_full_help() -> Manual { &[&monospace_roff("variance"), "Variance of coverage depths"], &[&monospace_roff("length"), "Length of each contig in base pairs"], &[&monospace_roff("count"), "Number of reads aligned to each contig. Note that supplementary alignments are not counted."], - &[&monospace_roff("metabat"), "(\"MetaBAT adjusted coverage\") Coverage as defined in Kang et al 2015 https://doi.org/10.7717/peerj.1165"], + &[&monospace_roff("metabat"), &format!("(\"MetaBAT adjusted coverage\") Coverage as defined in Kang et al 2015 https://doi.org/10.7717/peerj.1165. By default, reads with identity not strictly greater than 97% are excluded, matching MetaBAT's {}. A given {} uses CoverM's {} instead, so {} is not identical to the default. This method includes secondary alignments, and by default improper pair and supplementary alignments too; use {} or {} to exclude the latter two.", + &monospace_roff(">"), + &monospace_roff("--min-read-percent-identity"), + &monospace_roff(">="), + &monospace_roff("97"), + &monospace_roff("--proper-pairs-only"), + &monospace_roff("--exclude-supplementary"))], &[&monospace_roff("reads_per_base"), "Number of reads aligned divided by the length of the contig"], &[&monospace_roff("anir"), "Average BLAST-like identity of mapped reads"], &[&monospace_roff("rpkm"), "Reads mapped per kilobase of contig, per million mapped reads"], diff --git a/tests/test_cmdline.rs b/tests/test_cmdline.rs index 962f440..91f1d87 100644 --- a/tests/test_cmdline.rs +++ b/tests/test_cmdline.rs @@ -1956,6 +1956,52 @@ k141_2005182 225 1.9333 1.9333 0.0631", .unwrap(); } + #[test] + fn test_metabat_given_percent_identity_is_kept() { + // Reads below 97% identity, excluded by the metabat default, pass a + // given 90% threshold. + Assert::main_binary() + .with_args(&[ + "contig", + "-m", + "metabat", + "--min-read-percent-identity", + "90", + "-b", + "tests/data/k141_2005182.head11.bam", + ]) + .succeeds() + .stdout() + .contains( + "contigName contigLen totalAvgDepth k141_2005182.head11.bam k141_2005182.head11.bam-var +k141_2005182 225 6.04 6.04 1.12", + ) + .unwrap(); + } + + #[test] + fn test_metabat_given_97_uses_greater_or_equal() { + // A given 97 uses CoverM's >=, unlike the metabat default which + // approximates MetaBAT's > via 0.97001. Same BAM is 1.9333 by default. + Assert::main_binary() + .with_args(&[ + "contig", + "-m", + "metabat", + "--min-read-percent-identity", + "97", + "-b", + "tests/data/k141_2005182.head11.bam", + ]) + .succeeds() + .stdout() + .contains( + "contigName contigLen totalAvgDepth k141_2005182.head11.bam k141_2005182.head11.bam-var +k141_2005182 225 2.9333 2.9333 0.0631", + ) + .unwrap(); + } + #[test] fn test_deletions_count_towards_perc_id() { Assert::main_binary()