Mercurial > repos > iuc > remove_terminal_stop_codons

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/remove_terminal_stop_codons.py	Fri Dec 05 23:22:35 2025 +0000
@@ -0,0 +1,130 @@
+#!/usr/bin/env python3
+"""
+Remove terminal stop codons from coding sequences.
+
+Trim all terminal stop codons from sequences in a FASTA file, using a chosen
+NCBI genetic code (translation table). Leave non-stop terminal codons alone.
+If any INTERNAL, in-frame stop codon is found, exit with an error.
+
+This tool is designed as a preprocessing step for tools like cawlign and HyPhy
+that do not permit internal stop codons in their input sequences.
+
+Requires: Biopython
+"""
+
+import argparse
+import sys
+
+from Bio import SeqIO
+from Bio.Data import CodonTable
+from Bio.Seq import Seq
+
+
+def load_table(table_arg):
+    """Return a DNA codon table from an NCBI table id (int) or name (str)."""
+    if table_arg is None:
+        return CodonTable.unambiguous_dna_by_id[1]  # Standard
+    # try as integer id
+    try:
+        tid = int(table_arg)
+        return CodonTable.unambiguous_dna_by_id[tid]
+    except (ValueError, KeyError):
+        pass
+    # try as name
+    try:
+        return CodonTable.unambiguous_dna_by_name[table_arg]
+    except KeyError:
+        # Build a helpful hint list
+        valid_ids = sorted(CodonTable.unambiguous_dna_by_id.keys())
+        valid_names = sorted(CodonTable.unambiguous_dna_by_name.keys())
+        sys.stderr.write(
+            f"ERROR: Unknown genetic code '{table_arg}'.\n"
+            f"Try an NCBI table id (e.g., 1) or one of these names:\n"
+            f"  {', '.join(valid_names)}\n"
+            f"(Valid ids include: {', '.join(map(str, valid_ids))})\n"
+        )
+        sys.exit(2)
+
+
+def trim_terminal_stops_and_validate(record, stop_codons, check_internal=True):
+    """
+    Remove ALL trailing stop codons (0+ at the end).
+
+    If check_internal is True and any internal in-frame stop codon exists
+    (excluding the trailing block), exit with an error message.
+
+    Ignore a terminal codon that is not a stop codon.
+    """
+    # Work with DNA letters; treat any RNA U as T
+    seq_str = str(record.seq).upper().replace("U", "T")
+
+    # Count how many full codons sit at the end that are stops
+    idx = len(seq_str)
+    trailing_stops = 0
+    while idx >= 3:
+        codon = seq_str[idx - 3:idx]
+        if codon in stop_codons:
+            trailing_stops += 1
+            idx -= 3
+        else:
+            break
+
+    # Scan for INTERNAL stops: all complete codons up to (but not including)
+    # the trailing stop block (and ignoring any trailing partial codon).
+    if check_internal:
+        scan_end = (idx // 3) * 3  # only complete codons
+        for pos in range(0, scan_end, 3):
+            codon = seq_str[pos:pos + 3]
+            if codon in stop_codons:
+                sys.stderr.write(
+                    f"ERROR: Found an internal stop codon in sequence "
+                    f"'{record.id}' at position {pos}.\n"
+                    f"Tools like HyPhy and cawlign do not permit internal "
+                    f"stop codons. Please review your input sequences.\n"
+                )
+                sys.exit(2)
+
+    # Finally, remove the trailing stop codons (if any)
+    if trailing_stops > 0:
+        seq_str = seq_str[:idx]
+
+    # Leave sequences with non-stop terminal codons unchanged by design
+    return Seq(seq_str)
+
+
+def main():
+    ap = argparse.ArgumentParser(
+        description="Remove all terminal stop codons from a FASTA, using a "
+                    "chosen genetic code. Optionally fail if any internal "
+                    "in-frame stop codon is present."
+    )
+    ap.add_argument("-i", "--input", required=True, help="Input FASTA file")
+    ap.add_argument("-o", "--output", required=True, help="Output FASTA file")
+    ap.add_argument(
+        "-t", "--table",
+        help="NCBI translation table id (e.g., 1) or name "
+             "(e.g., 'Vertebrate Mitochondrial'). Default: 1 (Standard)."
+    )
+    ap.add_argument(
+        "--no-check-internal",
+        action="store_true",
+        help="Do not check for internal stop codons (only remove terminal)."
+    )
+    args = ap.parse_args()
+
+    table = load_table(args.table)
+    stop_codons = set(table.stop_codons)  # e.g., {'TAA','TAG','TGA'} for Standard
+
+    check_internal = not args.no_check_internal
+
+    records_out = []
+    for rec in SeqIO.parse(args.input, "fasta"):
+        new_seq = trim_terminal_stops_and_validate(rec, stop_codons, check_internal)
+        rec.seq = new_seq
+        records_out.append(rec)
+
+    SeqIO.write(records_out, args.output, "fasta")
+
+
+if __name__ == "__main__":
+    main()
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/remove_terminal_stop_codons.xml	Fri Dec 05 23:22:35 2025 +0000
@@ -0,0 +1,151 @@
+<tool id="remove_terminal_stop_codons" name="Remove terminal stop codons" version="@TOOL_VERSION@+galaxy@VERSION_SUFFIX@" profile="@PROFILE@">
+    <description>from coding sequences</description>
+    <macros>
+        <token name="@TOOL_VERSION@">1.0.0</token>
+        <token name="@VERSION_SUFFIX@">0</token>
+        <token name="@PROFILE@">25.0</token>
+    </macros>
+    <requirements>
+        <requirement type="package" version="1.84">biopython</requirement>
+    </requirements>
+    <required_files>
+        <include path="remove_terminal_stop_codons.py" />
+    </required_files>
+    <command detect_errors="exit_code"><![CDATA[
+        python '$__tool_directory__/remove_terminal_stop_codons.py'
+            -i '$input'
+            -o '$output'
+            #if str($genetic_code) != "1"
+                -t '$genetic_code'
+            #end if
+            $no_check_internal
+    ]]></command>
+    <inputs>
+        <param name="input" type="data" format="fasta" label="Input FASTA file"
+            help="FASTA file containing coding sequences (CDS) to process." />
+        <param name="genetic_code" type="select" label="Genetic code"
+            help="NCBI translation table to use for identifying stop codons. Different organisms use different genetic codes.">
+            <option value="1" selected="true">1 - Standard</option>
+            <option value="2">2 - Vertebrate Mitochondrial</option>
+            <option value="3">3 - Yeast Mitochondrial</option>
+            <option value="4">4 - Mold, Protozoan, Coelenterate Mitochondrial and Mycoplasma/Spiroplasma</option>
+            <option value="5">5 - Invertebrate Mitochondrial</option>
+            <option value="6">6 - Ciliate, Dasycladacean and Hexamita Nuclear</option>
+            <option value="9">9 - Echinoderm and Flatworm Mitochondrial</option>
+            <option value="10">10 - Euplotid Nuclear</option>
+            <option value="11">11 - Bacterial, Archaeal and Plant Plastid</option>
+            <option value="12">12 - Alternative Yeast Nuclear</option>
+            <option value="13">13 - Ascidian Mitochondrial</option>
+            <option value="14">14 - Alternative Flatworm Mitochondrial</option>
+            <option value="15">15 - Blepharisma Nuclear</option>
+            <option value="16">16 - Chlorophycean Mitochondrial</option>
+            <option value="21">21 - Trematode Mitochondrial</option>
+            <option value="22">22 - Scenedesmus obliquus Mitochondrial</option>
+            <option value="23">23 - Thraustochytrium Mitochondrial</option>
+            <option value="24">24 - Rhabdopleuridae Mitochondrial</option>
+            <option value="25">25 - Candidate Division SR1 and Gracilibacteria</option>
+            <option value="26">26 - Pachysolen tannophilus Nuclear</option>
+            <option value="27">27 - Karyorelict Nuclear</option>
+            <option value="28">28 - Condylostoma Nuclear</option>
+            <option value="29">29 - Mesodinium Nuclear</option>
+            <option value="30">30 - Peritrich Nuclear</option>
+            <option value="31">31 - Blastocrithidia Nuclear</option>
+            <option value="33">33 - Cephalodiscidae Mitochondrial UAA-Tyr</option>
+        </param>
+        <param argument="--no-check-internal" type="boolean" truevalue="--no-check-internal" falsevalue=""
+            checked="false" label="Skip internal stop codon check"
+            help="By default, the tool will fail if internal (in-frame) stop codons are found. Enable this to only remove terminal stops without checking for internal ones." />
+    </inputs>
+    <outputs>
+        <data name="output" format="fasta" label="${tool.name} on ${on_string}" />
+    </outputs>
+    <tests>
+        <!-- Test 1: Basic removal of terminal stop codon -->
+        <test expect_num_outputs="1">
+            <param name="input" value="with_terminal_stop.fasta" ftype="fasta" />
+            <param name="genetic_code" value="1" />
+            <output name="output" file="without_terminal_stop.fasta" ftype="fasta" />
+        </test>
+        <!-- Test 2: Sequence without terminal stop (should pass through unchanged) -->
+        <test expect_num_outputs="1">
+            <param name="input" value="no_terminal_stop.fasta" ftype="fasta" />
+            <param name="genetic_code" value="1" />
+            <output name="output" file="no_terminal_stop.fasta" ftype="fasta" />
+        </test>
+        <!-- Test 3: Internal stop codon should fail -->
+        <test expect_failure="true">
+            <param name="input" value="with_internal_stop.fasta" ftype="fasta" />
+            <param name="genetic_code" value="1" />
+        </test>
+        <!-- Test 4: Internal stop with skip check should pass -->
+        <test expect_num_outputs="1">
+            <param name="input" value="with_internal_stop.fasta" ftype="fasta" />
+            <param name="genetic_code" value="1" />
+            <param name="no_check_internal" value="true" />
+            <output name="output" file="with_internal_stop_output.fasta" ftype="fasta" />
+        </test>
+    </tests>
+    <help><![CDATA[
+**What it does**
+
+This tool removes terminal (trailing) stop codons from coding sequences in a FASTA file.
+It is designed as a **preprocessing step** for downstream tools like **cawlign** and **HyPhy**
+that do not permit stop codons in their input sequences.
+
+**Important**: By default, this tool will **fail with an error** if it detects any internal
+(in-frame) stop codons in your sequences. This is intentional, but can be disabled with the
+`--no-check-internal` option.
+
+----
+
+**Input**
+
+A FASTA file containing coding sequences (CDS). Sequences should be:
+
+- In the correct reading frame (starting at position 1 of a codon)
+- DNA sequences (RNA sequences with U will be converted to T)
+
+----
+
+**Output**
+
+A FASTA file with terminal stop codons removed. The output preserves:
+
+- Sequence identifiers and descriptions
+- Sequences that don't end with stop codons (passed through unchanged)
+- Partial codons at the end (not removed)
+
+----
+
+**Genetic Codes**
+
+Different organisms use different genetic codes (translation tables) which define
+which codons are stop codons:
+
+- **Standard (1)**: TAA, TAG, TGA - used by most organisms
+- **Vertebrate Mitochondrial (2)**: TAA, TAG, AGA, AGG - mitochondria of vertebrates
+- **Bacterial/Archaeal (11)**: TAA, TAG, TGA - bacteria and archaea
+
+Select the appropriate genetic code for your organism to ensure correct stop codon identification.
+
+----
+
+**Use Cases**
+
+1. **Before cawlign**: Remove terminal stops from sequences before codon-aware alignment
+2. **Before HyPhy**: Prepare sequences for selection analysis (HyPhy methods like BUSTED, FEL, MEME)
+3. **Quality control**: Identify sequences with internal stop codons that may need review
+    ]]></help>
+    <citations>
+        <citation type="bibtex">
+@misc{capheine2025,
+    author = {Callan, Danielle and Verdonk, Hannah and Kosakovsky Pond, Sergei L.},
+    title = {CAPHEINE: A Comprehensive Automated Pipeline Using HyPhy for Evolutionary Inference with Nextflow},
+    year = {2025},
+    publisher = {GitHub},
+    url = {https://github.com/veg/CAPHEINE},
+    note = {Terminal stop-codon removal logic in this Galaxy tool is adapted from the CAPHEINE pipeline.}
+}
+        </citation>
+    </citations>
+</tool>
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/test-data/no_terminal_stop.fasta	Fri Dec 05 23:22:35 2025 +0000
@@ -0,0 +1,4 @@
+>seq1 sequence without terminal stop
+ATGAAACCCGGGAAA
+>seq2 another sequence without terminal stop
+ATGCCCAAAGGGCCC
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/test-data/with_internal_stop.fasta	Fri Dec 05 23:22:35 2025 +0000
@@ -0,0 +1,2 @@
+>seq1 sequence with internal TAA stop codon
+ATGTAACCCGGGTAA
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/test-data/with_internal_stop_output.fasta	Fri Dec 05 23:22:35 2025 +0000
@@ -0,0 +1,2 @@
+>seq1 sequence with internal TAA stop codon
+ATGTAACCCGGG
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/test-data/with_terminal_stop.fasta	Fri Dec 05 23:22:35 2025 +0000
@@ -0,0 +1,6 @@
+>seq1 test sequence with terminal TAA stop
+ATGAAACCCGGGTAA
+>seq2 test sequence with terminal TAG stop
+ATGCCCAAAGGGCCCAAATAG
+>seq3 test sequence with terminal TGA stop
+ATGGGGTTTAAACCCGGGTGA
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/test-data/without_terminal_stop.fasta	Fri Dec 05 23:22:35 2025 +0000
@@ -0,0 +1,6 @@
+>seq1 test sequence with terminal TAA stop
+ATGAAACCCGGG
+>seq2 test sequence with terminal TAG stop
+ATGCCCAAAGGGCCCAAA
+>seq3 test sequence with terminal TGA stop
+ATGGGGTTTAAACCCGGG