Merge pull request #87 from gbouras13/issues

fix: unicycler --spades_options without --unicycler_options; canu pacbio-hifi (#83, #84)

Author: gbouras13

src/plassembler/__init__.py

@@ -38,6 +38,7 @@
 from plassembler.utils.plass_class import Assembly, Plass
 from plassembler.utils.qc import chopper, copy_sr_fastq_file, fastp, gzip_file
 from plassembler.utils.run_canu import (  # make_blastdb,; process_blast_output,; run_blast,
+    canu_read_type_and_error_rate,
     filter_entropy,
     filter_entropy_fastqs,
     run_canu,
@@ -1575,19 +1576,11 @@ def long(
         plasmidfastqs: Path = Path(outdir) / "plasmid_long.fastq"
         extract_long_fastqs_fast(samfile, plasmidfastqs, threads)
 
-        # to set error rate
-        canu_nano_or_pacbio = "nanopore"
-
-        if pacbio_model != "nothing":
-            if pacbio_model == "pacbio-hifi":
-                canu_nano_or_pacbio = "pacbio-hifi"
-                corrected_error_rate = 0.005
-            else:
-                canu_nano_or_pacbio = "pacbio"
-                corrected_error_rate = 0.045
-        else:
-            canu_nano_or_pacbio = "nanopore"
-            # corrected error rate default will be 0.12
+        # map the validated --pacbio_model to canu's read-type flag, error rate,
+        # and whether read correction should be skipped (HiFi reads)
+        canu_nano_or_pacbio, corrected_error_rate, skip_canu_correct = (
+            canu_read_type_and_error_rate(pacbio_model, corrected_error_rate)
+        )
 
         if canu_flag is True:
             assembler = "canu"
@@ -1641,28 +1634,38 @@ def long(
             )
             filter_entropy_fastqs(plasmidfastqs, entropy_filtered_fastq)
 
-            logger.info("Correcting reads with canu prior to running Unicycler.")
-
-            try:
-                run_canu_correct(
-                    threads,
-                    logdir,
-                    entropy_filtered_fastq,
-                    canu_output_dir,
-                    canu_nano_or_pacbio,
-                    total_flye_plasmid_length,
-                    corrected_error_rate,
-                    coverage,
-                )
-                # convert the corrected .fasta.gz from Canu to fastq
-                canu_reads: Path = (
-                    Path(canu_output_dir) / "canu.correctedReads.fasta.gz"
+            if skip_canu_correct:
+                # PacBio HiFi reads are already high-accuracy; canu -correct
+                # rejects them ("Cannot correct already corrected reads"), so
+                # skip correction and assemble the entropy-filtered reads (#84).
+                logger.info(
+                    "PacBio HiFi reads detected - skipping canu read correction."
                 )
-                corrected_fastqs: Path = Path(outdir) / "corrected_plasmid_long.fastq"
-                corrected_fasta_to_fastq(canu_reads, corrected_fastqs)
-            except Exception:
-                logger.warning("Advancing with uncorrected reads")
                 corrected_fastqs = entropy_filtered_fastq
+            else:
+                logger.info("Correcting reads with canu prior to running Unicycler.")
+                try:
+                    run_canu_correct(
+                        threads,
+                        logdir,
+                        entropy_filtered_fastq,
+                        canu_output_dir,
+                        canu_nano_or_pacbio,
+                        total_flye_plasmid_length,
+                        corrected_error_rate,
+                        coverage,
+                    )
+                    # convert the corrected .fasta.gz from Canu to fastq
+                    canu_reads: Path = (
+                        Path(canu_output_dir) / "canu.correctedReads.fasta.gz"
+                    )
+                    corrected_fastqs: Path = (
+                        Path(outdir) / "corrected_plasmid_long.fastq"
+                    )
+                    corrected_fasta_to_fastq(canu_reads, corrected_fastqs)
+                except Exception:
+                    logger.warning("Advancing with uncorrected reads")
+                    corrected_fastqs = entropy_filtered_fastq
 
             # remove canu directory
             remove_directory(canu_output_dir)

src/plassembler/utils/run_canu.py

@@ -9,6 +9,29 @@
 from plassembler.utils.external_tools import ExternalTool
 
 
+def canu_read_type_and_error_rate(pacbio_model, default_error_rate):
+    """Map a validated ``--pacbio_model`` value to canu's read-type flag, the
+    corrected error rate, and whether canu read correction should be skipped.
+
+    ``pacbio_model`` is the output of ``validate_pacbio_model`` (so one of
+    ``"--pacbio-raw"``, ``"--pacbio-corr"``, ``"--pacbio-hifi"``) or
+    ``"nothing"``.
+
+    PacBio HiFi reads are already high-accuracy: canu's read type is
+    ``-pacbio-hifi`` and its correction step is skipped (running ``canu -correct``
+    on HiFi reads errors with "Cannot correct already corrected reads" - issue
+    #84). Previously the check compared against ``"pacbio-hifi"`` (without the
+    ``--`` prefix), so HiFi reads were wrongly assembled as regular ``-pacbio``.
+
+    :return: ``(canu_read_type, corrected_error_rate, skip_correction)``
+    """
+    if pacbio_model == "--pacbio-hifi":
+        return "pacbio-hifi", 0.005, True
+    if pacbio_model != "nothing":
+        return "pacbio", 0.045, False
+    return "nanopore", default_error_rate, False
+
+
 def run_canu_correct(
     threads: Path,
     logdir: Path,

src/plassembler/utils/run_unicycler.py

@@ -6,6 +6,22 @@
 from plassembler.utils.external_tools import ExternalTool
 
 
+def build_extra_unicycler_options(unicycler_options, spades_options):
+    """Build the extra-options suffix for a unicycler command.
+
+    Each of ``unicycler_options`` / ``spades_options`` may be ``None``
+    independently. Regression (issue #83): when only ``spades_options`` was set,
+    a ``None`` ``unicycler_options`` was interpolated into the command as the
+    literal string ``"None"``, which Unicycler rejected.
+    """
+    parts = []
+    if unicycler_options is not None:
+        parts.append(str(unicycler_options))
+    if spades_options is not None:
+        parts.append(f'--spades_options "{spades_options}"')
+    return " ".join(parts)
+
+
 def run_unicycler(
     threads: int,
     logdir: Path,
@@ -28,34 +44,18 @@ def run_unicycler(
     :return:
     """
 
-    if unicycler_options is None and spades_options is None:
-        unicycler = ExternalTool(
-            tool="unicycler",
-            input="",
-            output="",
-            params=f" -1 {short_one} -2 {short_two} -l {longreads} -t {threads} -o {unicycler_output_dir}",
-            logdir=logdir,
-            outfile="",
-        )
-    else:
-        if spades_options is None:
-            unicycler = ExternalTool(
-                tool="unicycler",
-                input="",
-                output="",
-                params=f" -1 {short_one} -2 {short_two} -l {longreads} -t {threads} -o {unicycler_output_dir} {unicycler_options}",
-                logdir=logdir,
-                outfile="",
-            )
-        else:
-            unicycler = ExternalTool(
-                tool="unicycler",
-                input="",
-                output="",
-                params=f' -1 {short_one} -2 {short_two} -l {longreads} -t {threads} -o {unicycler_output_dir} {unicycler_options} --spades_options "{spades_options}" ',
-                logdir=logdir,
-                outfile="",
-            )
+    extra = build_extra_unicycler_options(unicycler_options, spades_options)
+    unicycler = ExternalTool(
+        tool="unicycler",
+        input="",
+        output="",
+        params=(
+            f" -1 {short_one} -2 {short_two} -l {longreads} -t {threads} "
+            f"-o {unicycler_output_dir} {extra}"
+        ),
+        logdir=logdir,
+        outfile="",
+    )
 
     ExternalTool.run_tool(unicycler, to_stdout=False)
 
@@ -80,34 +80,18 @@ def run_unicycler_long(
     :return:
     """
 
-    if unicycler_options is None and spades_options is None:
-        unicycler_long = ExternalTool(
-            tool="unicycler",
-            input="",
-            output="",
-            params=f" -s {corrected_longreads} -l {entropy_filtered_longreads} -t {threads} -o {unicycler_output_dir}",
-            logdir=logdir,
-            outfile="",
-        )
-    else:
-        if spades_options is None:
-            unicycler_long = ExternalTool(
-                tool="unicycler",
-                input="",
-                output="",
-                params=f" -s {corrected_longreads} -l {entropy_filtered_longreads} -t {threads} -o {unicycler_output_dir} {unicycler_options}",
-                logdir=logdir,
-                outfile="",
-            )
-        else:
-            unicycler_long = ExternalTool(
-                tool="unicycler",
-                input="",
-                output="",
-                params=f' -s {corrected_longreads} -l {entropy_filtered_longreads} -t {threads} -o {unicycler_output_dir} {unicycler_options} --spades_options "{spades_options}"',
-                logdir=logdir,
-                outfile="",
-            )
+    extra = build_extra_unicycler_options(unicycler_options, spades_options)
+    unicycler_long = ExternalTool(
+        tool="unicycler",
+        input="",
+        output="",
+        params=(
+            f" -s {corrected_longreads} -l {entropy_filtered_longreads} "
+            f"-t {threads} -o {unicycler_output_dir} {extra}"
+        ),
+        logdir=logdir,
+        outfile="",
+    )
 
     ExternalTool.run_tool(unicycler_long, to_stdout=False)
 

tests/test_issues.py

@@ -0,0 +1,59 @@
+"""Regression tests for reported GitHub issues."""
+
+import pytest
+
+from src.plassembler.utils.run_canu import canu_read_type_and_error_rate
+from src.plassembler.utils.run_unicycler import build_extra_unicycler_options
+
+# ---------------------------------------------------------------------------
+# issue #83: Unicycler crashes with --spades_options but no --unicycler_options
+# (a None unicycler_options was interpolated into the command as "None")
+# ---------------------------------------------------------------------------
+
+
+def test_unicycler_options_both_none():
+    assert build_extra_unicycler_options(None, None) == ""
+
+
+def test_unicycler_options_spades_only_has_no_none():
+    out = build_extra_unicycler_options(None, "--memory 2000")
+    assert out == '--spades_options "--memory 2000"'
+    assert "None" not in out
+
+
+def test_unicycler_options_unicycler_only():
+    assert build_extra_unicycler_options("--mode bold", None) == "--mode bold"
+
+
+def test_unicycler_options_both():
+    out = build_extra_unicycler_options("--mode bold", "--memory 2000")
+    assert out == '--mode bold --spades_options "--memory 2000"'
+
+
+# ---------------------------------------------------------------------------
+# issue #84: `plassembler long --pacbio_model pacbio-hifi` passed the wrong
+# read type to canu (-pacbio instead of -pacbio-hifi) and tried to correct
+# already-corrected HiFi reads
+# ---------------------------------------------------------------------------
+
+
+def test_canu_read_type_hifi_skips_correction():
+    read_type, err, skip = canu_read_type_and_error_rate("--pacbio-hifi", 0.12)
+    assert read_type == "pacbio-hifi"
+    assert err == 0.005
+    assert skip is True
+
+
+@pytest.mark.parametrize("model", ["--pacbio-raw", "--pacbio-corr"])
+def test_canu_read_type_pacbio(model):
+    read_type, err, skip = canu_read_type_and_error_rate(model, 0.12)
+    assert read_type == "pacbio"
+    assert err == 0.045
+    assert skip is False
+
+
+def test_canu_read_type_nanopore_passthrough():
+    read_type, err, skip = canu_read_type_and_error_rate("nothing", 0.12)
+    assert read_type == "nanopore"
+    assert err == 0.12  # keeps the CLI --corrected_error_rate
+    assert skip is False