Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
69 changes: 36 additions & 33 deletions src/plassembler/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -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,
Expand Down Expand Up @@ -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"
Expand Down Expand Up @@ -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)
Expand Down
23 changes: 23 additions & 0 deletions src/plassembler/utils/run_canu.py
Original file line number Diff line number Diff line change
Expand Up @@ -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,
Expand Down
96 changes: 40 additions & 56 deletions src/plassembler/utils/run_unicycler.py
Original file line number Diff line number Diff line change
Expand Up @@ -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,
Expand All @@ -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)

Expand All @@ -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)

Expand Down
59 changes: 59 additions & 0 deletions tests/test_issues.py
Original file line number Diff line number Diff line change
@@ -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
Loading