From b1d6ce3d35026f0ac1f9abd0c1de3b594e78519c Mon Sep 17 00:00:00 2001 From: Robert Turnbull Date: Tue, 14 Oct 2025 12:11:26 +1100 Subject: [PATCH 01/11] :zap: allowing for custom training data --- terrier/apps.py | 38 ++++++++++--- terrier/repeatmasker.py | 83 ++++++++++++++++++--------- tests/test-data/custom.fna | 4 ++ tests/test-data/repbase.ref | 6 ++ tests/test-data/simple.ref | 4 ++ tests/test_repeatmasker.py | 109 ++++++++++++++++++++++++++++++++++++ 6 files changed, 207 insertions(+), 37 deletions(-) create mode 100644 tests/test-data/custom.fna create mode 100644 tests/test-data/repbase.ref create mode 100644 tests/test-data/simple.ref create mode 100644 tests/test_repeatmasker.py diff --git a/terrier/apps.py b/terrier/apps.py index 55c605d..948ebe9 100644 --- a/terrier/apps.py +++ b/terrier/apps.py @@ -233,7 +233,7 @@ def create_repeatmasker_seqtree(self, output:Path, repbase:Path, label_smoothing @ta.tool def preprocess( self, - repbase:Path=ta.Param(..., help="The path to the RepBase fasta directory."), + input:list[Path]=ta.Param(..., help="The path to a FASTA file, multiple FASTA files or a directory of FASTA files (e.g. the RepBase FASTA directory)"), seqbank:Path=ta.Param(..., help="The path to save the new SeqBank file."), seqtree:Path=ta.Param(..., help="The path to save the new SeqTree file."), label_smoothing:float=0.0, @@ -243,23 +243,43 @@ def preprocess( from seqbank import SeqBank from .repeatmasker import create_repeatmasker_seqtree + ALLOWABLE_EXTENSIONS = [ + "fasta", + "fa", + "fna", + "ref", + ] + seqbank = SeqBank(path=seqbank, write=True) - assert repbase is not None - repbase = Path(repbase) - assert repbase.exists() + + # File all the FASTA files + fasta_paths = [] + for input_path in input: + input_path = Path(input_path) + if input_path.is_dir(): + for extension in ALLOWABLE_EXTENSIONS: + fasta_paths += list(input_path.glob(f'*.{extension}')) + elif input_path.is_file() and input_path.exists(): + fasta_paths.append(input_path) + assert len(fasta_paths), ( + f"No FASTA files found. " + f"Make sure you provide at least one file or a directory of files with an expected FASTA extension:" + f"({', '.join(ALLOWABLE_EXTENSIONS)})" + ) # Create the seqbank from the FASTA files with .ref extension - files = list(repbase.glob('*.ref')) - seqbank.add_files(files, format="fasta") + seqbank.add_files(fasta_paths, format="fasta") # Create the seqtree - return create_repeatmasker_seqtree( - output=seqtree, - repbase=repbase, + seqtree_path = Path(seqtree) + seqtree = create_repeatmasker_seqtree( + paths=fasta_paths, label_smoothing=label_smoothing, gamma=gamma, partitions=partitions, ) + seqtree.save(seqtree_path) + seqtree.classification_tree.render(print=1) @ta.tool def evaluate( diff --git a/terrier/repeatmasker.py b/terrier/repeatmasker.py index c7e6108..0638319 100644 --- a/terrier/repeatmasker.py +++ b/terrier/repeatmasker.py @@ -6,7 +6,32 @@ from collections import Counter -def create_repeatmasker_seqtree(output:Path, repbase:Path, label_smoothing:float=0.0, gamma:float=0.0, partitions:int=5): +def get_verbatim_classification(path:Path, record) -> str: + accession = record.id + + description = record.description + + if "#" in description: + description = description[ description.find("#")+1 : ] + + components = description.split("\t") + if len(components) == 3: + return components[1] + + if path.name == "simple.ref": + return "Simple Repeat" + elif accession.startswith("SINE_"): + return "SINE" + else: + return description + + +def create_repeatmasker_seqtree( + fasta_paths:list[Path], + label_smoothing:float=0.0, + gamma:float=0.0, + partitions:int=5, +) -> SeqTree: with open(Path(__file__).parent/"data/repbase-to-repeatmasker.toml", "r") as f: mapping = toml.load(f) @@ -21,45 +46,49 @@ def create_repeatmasker_seqtree(output:Path, repbase:Path, label_smoothing:float # Read files count = 0 - for file in repbase.glob('*.ref'): + for file in fasta_paths: with open(file) as f: for record in SeqIO.parse(f, "fasta"): partition = count % partitions accession = record.id - components = record.description.split("\t") - if len(components) != 3: - if file.name == "simple.ref": - classification = "Simple Repeat" - elif accession.startswith("SINE_"): - classification = "SINE" - else: - continue - else: - classification = components[1] - - if classification not in mapping: + classification = get_verbatim_classification(file, record) + + if classification in mapping: + mapped_counter.update([classification]) + classification = mapping[classification] + + if classification not in mapping.values(): not_mapped_counter.update([classification]) continue - mapped_counter.update([classification]) - - repeat_name = mapping[classification] - if repeat_name == "Unknown": + if classification == "Unknown": continue - if repeat_name not in classification_nodes: - components = repeat_name.split("/") + if classification not in classification_nodes: + components = classification.split("/") repeat_type = components[0] repeat_subtype = components[1] if len(components) > 1 else "" if repeat_type not in classification_nodes: - classification_nodes[repeat_type] = SoftmaxNode(repeat_type, parent=classification_tree, label_smoothing=label_smoothing, gamma=gamma, repeat_masker_name=repeat_type) + classification_nodes[repeat_type] = SoftmaxNode( + repeat_type, + parent=classification_tree, + label_smoothing=label_smoothing, + gamma=gamma, + repeat_masker_name=repeat_type, + ) repeat_type_node = classification_nodes[repeat_type] if repeat_subtype: - classification_nodes[repeat_name] = SoftmaxNode(repeat_subtype, parent=repeat_type_node, label_smoothing=label_smoothing, gamma=gamma, repeat_masker_name=repeat_name) - - node = classification_nodes[repeat_name] + classification_nodes[classification] = SoftmaxNode( + repeat_subtype, + parent=repeat_type_node, + label_smoothing=label_smoothing, + gamma=gamma, + repeat_masker_name=classification, + ) + + node = classification_nodes[classification] try: seqtree.add(accession, node, partition) @@ -68,13 +97,11 @@ def create_repeatmasker_seqtree(output:Path, repbase:Path, label_smoothing:float count += 1 - print("repbase,count,mapped,repeat_masker") + print("provided,count,mapped,repeat_masker") for classification,count in mapped_counter.most_common(): print(classification,count,1, mapping[classification], sep=",") for classification,count in not_mapped_counter.most_common(): print(classification,count,0, "", sep=",") - seqtree.save(output) - seqtree.classification_tree.render(print=1) - + return seqtree diff --git a/tests/test-data/custom.fna b/tests/test-data/custom.fna new file mode 100644 index 0000000..fe6abc5 --- /dev/null +++ b/tests/test-data/custom.fna @@ -0,0 +1,4 @@ +>I123#DNA/Academ +ACGT +>XXadsfadsf#LTR/Caulimovirus +AGCT \ No newline at end of file diff --git a/tests/test-data/repbase.ref b/tests/test-data/repbase.ref new file mode 100644 index 0000000..6cb5718 --- /dev/null +++ b/tests/test-data/repbase.ref @@ -0,0 +1,6 @@ +>IS905 DNA transposon Lactococcus lactis +ACGT +>BAGGINS1 Loa Drosophila melanogaster +ACGT +>SINE_DFSFDFS +ACGT \ No newline at end of file diff --git a/tests/test-data/simple.ref b/tests/test-data/simple.ref new file mode 100644 index 0000000..c243294 --- /dev/null +++ b/tests/test-data/simple.ref @@ -0,0 +1,4 @@ +>IJFDKSFF +ACGT +>GOOGLE +ACGT \ No newline at end of file diff --git a/tests/test_repeatmasker.py b/tests/test_repeatmasker.py new file mode 100644 index 0000000..82449fe --- /dev/null +++ b/tests/test_repeatmasker.py @@ -0,0 +1,109 @@ +from pathlib import Path +from Bio import SeqIO +from terrier.repeatmasker import get_verbatim_classification, create_repeatmasker_seqtree + + +TEST_DATA_DIR = Path(__file__).parent / "test-data" + +def test_get_verbatim_classification_repbase(): + fasta = TEST_DATA_DIR / "repbase.ref" + expected = [ + "DNA transposon", + "Loa", + "SINE", + ] + count = sum(1 for line in open(fasta) if line.startswith(">")) + assert count == len(expected) + with open(fasta) as f: + for record, exp in zip(SeqIO.parse(f, "fasta"), expected): + classification = get_verbatim_classification(fasta, record) + assert classification == exp + + +def test_get_verbatim_classification_simple(): + fasta = TEST_DATA_DIR / "simple.ref" + expected = [ + "Simple Repeat", + "Simple Repeat", + ] + count = sum(1 for line in open(fasta) if line.startswith(">")) + assert count == len(expected) + with open(fasta) as f: + for record, exp in zip(SeqIO.parse(f, "fasta"), expected): + classification = get_verbatim_classification(fasta, record) + assert classification == exp + + +def test_get_verbatim_classification_custom(): + expected = [ + "DNA/Academ", + "LTR/Caulimovirus", + ] + fasta = TEST_DATA_DIR / "custom.fna" + + count = sum(1 for line in open(fasta) if line.startswith(">")) + assert count == len(expected) + with open(fasta) as f: + for record, exp in zip(SeqIO.parse(f, "fasta"), expected): + classification = get_verbatim_classification(fasta, record) + assert classification == exp + + +def test_create_repeatmasker_seqtree_repbase(): + fasta = TEST_DATA_DIR / "repbase.ref" + seqtree = create_repeatmasker_seqtree([fasta]) + assert seqtree.classification_tree is not None + assert seqtree is not None + assert len(seqtree) == 3 + assert {"IS905", "BAGGINS1", "SINE_DFSFDFS"} == set(seqtree.keys()) + assert seqtree.classification_tree.render_equal( + """ + root + ├── DNA + ├── LINE + │ └── R1 + └── SINE + """ + ) + +def test_create_repeatmasker_seqtree_repbase_simple(): + seqtree = create_repeatmasker_seqtree([TEST_DATA_DIR / "repbase.ref", TEST_DATA_DIR / "simple.ref"]) + assert seqtree.classification_tree is not None + assert seqtree is not None + assert len(seqtree) == 5 + assert {"IS905", "BAGGINS1", "SINE_DFSFDFS", "IJFDKSFF", "GOOGLE"} == set(seqtree.keys()) + assert seqtree.classification_tree.render_equal( + """ + root + ├── DNA + ├── LINE + │ └── R1 + ├── SINE + └── Satellite + """ + ) + +def test_create_repeatmasker_seqtree_custom(): + seqtree = create_repeatmasker_seqtree([ + TEST_DATA_DIR / "repbase.ref", + TEST_DATA_DIR / "simple.ref", + TEST_DATA_DIR / "custom.fna", + ]) + assert seqtree.classification_tree is not None + assert seqtree is not None + assert len(seqtree) == 7 + assert {"IS905", "BAGGINS1", "SINE_DFSFDFS", "IJFDKSFF", "GOOGLE", "I123#DNA/Academ", "XXadsfadsf#LTR/Caulimovirus"} == set(seqtree.keys()) + assert seqtree.classification_tree.render_equal( + """ + root + ├── DNA + │ └── Academ + ├── LINE + │ └── R1 + ├── SINE + ├── Satellite + └── LTR + └── Caulimovirus + """ + ) + assert seqtree["I123#DNA/Academ"].node.name == "Academ" From a8e754ac25f9d0da04c589cfdda4b5e017a125ce Mon Sep 17 00:00:00 2001 From: Robert Turnbull Date: Tue, 14 Oct 2025 14:10:02 +1100 Subject: [PATCH 02/11] :bug: fixing way of building the seqtree --- docs/preprocessing.rst | 37 ++++++++++++++++-- docs/training.rst | 13 +++++++ terrier/apps.py | 85 +++++++++++++++++++---------------------- terrier/repeatmasker.py | 2 +- 4 files changed, 87 insertions(+), 50 deletions(-) diff --git a/docs/preprocessing.rst b/docs/preprocessing.rst index 44870d5..8f71575 100644 --- a/docs/preprocessing.rst +++ b/docs/preprocessing.rst @@ -28,7 +28,7 @@ These two files can be generated from the Repbase database using the ``terrier-t .. code-block:: bash - terrier-tools preprocess --repbase $REPBASE_DIR --seqbank $REPBASE_DIR/Repbase-seqbank.sb --seqtree $REPBASE_DIR/Repbase-seqtree.st + terrier-tools preprocess --input $REPBASE_DIR --seqbank $REPBASE_DIR/Repbase-seqbank.sb --seqtree $REPBASE_DIR/Repbase-seqtree.st This will create a SeqBank file called ``Repbase-seqbank.sb`` and a SeqTree files called ``Repbase-seqtree.st`` and place them the ``$REPBASE_DIR``. @@ -46,13 +46,13 @@ Now you are ready to train Terrier using the SeqBank and SeqTree files you have Optional: Display the SeqTree ------------------------------ -You can list the number of accessions for each node in the SeqTree file with this command: +You can list the number of sequences for each node in the SeqTree file with this command: .. code-block:: bash seqtree render $REPBASE_DIR/Repbase-seqtree.st --print --count -That will output a tree with the number of accessions like this: +That will output a tree with the number of sequences like this: .. code-block:: text @@ -124,4 +124,33 @@ This will create an HTML file with the Sunburst chart of the SeqTree like this: You can open the HTML file in a browser to view the chart. -You can also output the SeqTree with a .png, .svg, or .pdf extension by changing the extension of the output file. \ No newline at end of file +You can also output the SeqTree with a .png, .svg, or .pdf extension by changing the extension of the output file. + +Custom Datasets +---------------- + +You can create a custom repeat library in FASTA format, with the classification of each sequence like this: + +.. code-block:: text + + >SeqID#DNA/Academ + ACTGACTGACTG... + +Or with the classification separated with a tab character like this: + +.. code-block:: text + + >SeqID LTR/Caulimovirus + ACTGACTGACTG... + +Then preprocess like this: + +.. code-block:: bash + + terrier-tools preprocess --input custom.fasta --seqbank custom-seqbank.sb --seqtree custom-seqtree.st + +You can include Repbase with your custom dataset like this: + +.. code-block:: bash + + terrier-tools preprocess --input $REPBASE_DIR --input custom.fasta --seqbank combined-seqbank.sb --seqtree combined-seqtree.st \ No newline at end of file diff --git a/docs/training.rst b/docs/training.rst index 3a7bc82..59d64d9 100644 --- a/docs/training.rst +++ b/docs/training.rst @@ -17,6 +17,19 @@ To use the same hyperparameters as in the main release of Terrier, you can run t --seqtree $SEQTREE \ --seqbank $SEQBANK +If you want to train using the pretrained Terrier model weights as a starting point, you can add the ``--pretrained`` flag: + +.. code-block:: bash + + SEQBANK=$REPBASE_DIR/Repbase-seqbank.sb + SEQTREE=$REPBASE_DIR/Repbase-seqtree.st + terrier-tools train \ + --seqtree $SEQTREE \ + --seqbank $SEQBANK \ + --pretrained default + +You can replace the word ``default`` with a path to a checkpoint file if you have one or to a URL to a checkpoint file. + You can see other command-line options by running: .. code-block:: bash diff --git a/terrier/apps.py b/terrier/apps.py index 948ebe9..238f6fd 100644 --- a/terrier/apps.py +++ b/terrier/apps.py @@ -219,67 +219,62 @@ def get_prediction_probability(row): def checkpoint(self, checkpoint:Path=None) -> str: return checkpoint or "https://github.com/rbturnbull/terrier/releases/download/v0.2.0/terrier-0.2.0.ckpt" - @ta.tool - def create_repeatmasker_seqtree(self, output:Path, repbase:Path, label_smoothing:float=0.0, gamma:float=0.0, partitions:int=5): - from .repeatmasker import create_repeatmasker_seqtree - return create_repeatmasker_seqtree( - output=output, - repbase=repbase, - label_smoothing=label_smoothing, - gamma=gamma, - partitions=partitions, - ) - - @ta.tool - def preprocess( - self, - input:list[Path]=ta.Param(..., help="The path to a FASTA file, multiple FASTA files or a directory of FASTA files (e.g. the RepBase FASTA directory)"), - seqbank:Path=ta.Param(..., help="The path to save the new SeqBank file."), - seqtree:Path=ta.Param(..., help="The path to save the new SeqTree file."), - label_smoothing:float=0.0, - gamma:float=0.0, - partitions:int=5, - ): - from seqbank import SeqBank - from .repeatmasker import create_repeatmasker_seqtree - + def find_fasta_paths(self, files:list[Path]) -> list[Path]: ALLOWABLE_EXTENSIONS = [ "fasta", "fa", "fna", "ref", ] - - seqbank = SeqBank(path=seqbank, write=True) - - # File all the FASTA files fasta_paths = [] - for input_path in input: - input_path = Path(input_path) - if input_path.is_dir(): + for file in files: + file = Path(file) + if file.is_dir(): for extension in ALLOWABLE_EXTENSIONS: - fasta_paths += list(input_path.glob(f'*.{extension}')) - elif input_path.is_file() and input_path.exists(): - fasta_paths.append(input_path) + fasta_paths += list(file.glob(f'*.{extension}')) + elif file.is_file() and file.exists(): + fasta_paths.append(file) assert len(fasta_paths), ( f"No FASTA files found. " f"Make sure you provide at least one file or a directory of files with an expected FASTA extension:" f"({', '.join(ALLOWABLE_EXTENSIONS)})" ) + fasta_paths = [self.process_location(file) for file in fasta_paths] + return fasta_paths - # Create the seqbank from the FASTA files with .ref extension - seqbank.add_files(fasta_paths, format="fasta") + @ta.tool("create_repeatmasker_seqtree") + def preprocess( + self, + input:list[Path]=ta.Param(..., help="The path to a FASTA file, multiple FASTA files or a directory of FASTA files (e.g. the RepBase FASTA directory)"), + seqbank:Path=ta.Param(None, help="The path to save the new SeqBank file."), + seqtree:Path=ta.Param(None, help="The path to save the new SeqTree file."), + label_smoothing:float=0.0, + gamma:float=0.0, + partitions:int=5, + ): + from seqbank import SeqBank + from .repeatmasker import create_repeatmasker_seqtree + + assert seqbank or seqtree, "You must provide either a --seqbank or --seqtree output path (usually both)." + + fasta_paths = self.find_fasta_paths(input) + + # Create the seqbank from the FASTA files + if seqbank: + seqbank = SeqBank(path=seqbank, write=True) + seqbank.add_files(fasta_paths, format="fasta") # Create the seqtree - seqtree_path = Path(seqtree) - seqtree = create_repeatmasker_seqtree( - paths=fasta_paths, - label_smoothing=label_smoothing, - gamma=gamma, - partitions=partitions, - ) - seqtree.save(seqtree_path) - seqtree.classification_tree.render(print=1) + if seqtree: + seqtree_path = Path(seqtree) + seqtree = create_repeatmasker_seqtree( + fasta_paths=fasta_paths, + label_smoothing=label_smoothing, + gamma=gamma, + partitions=partitions, + ) + seqtree.save(seqtree_path) + seqtree.classification_tree.render(print=1) @ta.tool def evaluate( diff --git a/terrier/repeatmasker.py b/terrier/repeatmasker.py index 0638319..234a69e 100644 --- a/terrier/repeatmasker.py +++ b/terrier/repeatmasker.py @@ -15,7 +15,7 @@ def get_verbatim_classification(path:Path, record) -> str: description = description[ description.find("#")+1 : ] components = description.split("\t") - if len(components) == 3: + if len(components) >= 2: return components[1] if path.name == "simple.ref": From adb0feae880bd82813e696183dbb73898c557d47 Mon Sep 17 00:00:00 2001 From: Robert Turnbull Date: Tue, 14 Oct 2025 14:11:07 +1100 Subject: [PATCH 03/11] :zap: removing reference to extra tool --- terrier/apps.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/terrier/apps.py b/terrier/apps.py index 238f6fd..8d4a7b9 100644 --- a/terrier/apps.py +++ b/terrier/apps.py @@ -242,7 +242,7 @@ def find_fasta_paths(self, files:list[Path]) -> list[Path]: fasta_paths = [self.process_location(file) for file in fasta_paths] return fasta_paths - @ta.tool("create_repeatmasker_seqtree") + @ta.tool def preprocess( self, input:list[Path]=ta.Param(..., help="The path to a FASTA file, multiple FASTA files or a directory of FASTA files (e.g. the RepBase FASTA directory)"), From 92e720c0bb523a4866999ab67f70427f424b279c Mon Sep 17 00:00:00 2001 From: Robert Turnbull Date: Tue, 14 Oct 2025 14:15:18 +1100 Subject: [PATCH 04/11] :zap: showing the counts in the tree --- terrier/apps.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/terrier/apps.py b/terrier/apps.py index 8d4a7b9..2e32140 100644 --- a/terrier/apps.py +++ b/terrier/apps.py @@ -274,7 +274,7 @@ def preprocess( partitions=partitions, ) seqtree.save(seqtree_path) - seqtree.classification_tree.render(print=1) + seqtree.classification_tree.render(print=1, count=True) @ta.tool def evaluate( From 5f21ad89fa7475d13ebbf5b6dbec3d4228d496d0 Mon Sep 17 00:00:00 2001 From: Robert Turnbull Date: Tue, 14 Oct 2025 14:16:08 +1100 Subject: [PATCH 05/11] :bug: fixing issue rendering tree --- terrier/apps.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/terrier/apps.py b/terrier/apps.py index 2e32140..46e1f3f 100644 --- a/terrier/apps.py +++ b/terrier/apps.py @@ -274,7 +274,7 @@ def preprocess( partitions=partitions, ) seqtree.save(seqtree_path) - seqtree.classification_tree.render(print=1, count=True) + seqtree.render(print=1, count=True) @ta.tool def evaluate( From 72e1fc82a8b825fad4ded789246135016aa3a011 Mon Sep 17 00:00:00 2001 From: rturnbull Date: Tue, 14 Oct 2025 14:29:06 +1100 Subject: [PATCH 06/11] :zap: allowing URLs for preprocessing --- terrier/apps.py | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/terrier/apps.py b/terrier/apps.py index 46e1f3f..4264506 100644 --- a/terrier/apps.py +++ b/terrier/apps.py @@ -219,16 +219,18 @@ def get_prediction_probability(row): def checkpoint(self, checkpoint:Path=None) -> str: return checkpoint or "https://github.com/rbturnbull/terrier/releases/download/v0.2.0/terrier-0.2.0.ckpt" - def find_fasta_paths(self, files:list[Path]) -> list[Path]: + def find_fasta_paths(self, files:list[str]) -> list[Path]: ALLOWABLE_EXTENSIONS = [ "fasta", "fa", "fna", + "fas", "ref", ] fasta_paths = [] for file in files: - file = Path(file) + file = self.process_location(file) + if file.is_dir(): for extension in ALLOWABLE_EXTENSIONS: fasta_paths += list(file.glob(f'*.{extension}')) @@ -239,13 +241,12 @@ def find_fasta_paths(self, files:list[Path]) -> list[Path]: f"Make sure you provide at least one file or a directory of files with an expected FASTA extension:" f"({', '.join(ALLOWABLE_EXTENSIONS)})" ) - fasta_paths = [self.process_location(file) for file in fasta_paths] return fasta_paths @ta.tool def preprocess( self, - input:list[Path]=ta.Param(..., help="The path to a FASTA file, multiple FASTA files or a directory of FASTA files (e.g. the RepBase FASTA directory)"), + input:list[str]=ta.Param(..., help="The path to a FASTA file, URL to a FASTA file, multiple FASTA files or a directory of FASTA files (e.g. the RepBase FASTA directory)"), seqbank:Path=ta.Param(None, help="The path to save the new SeqBank file."), seqtree:Path=ta.Param(None, help="The path to save the new SeqTree file."), label_smoothing:float=0.0, From 0f41bca574becd55d227104e1e71cc99359d838f Mon Sep 17 00:00:00 2001 From: rturnbull Date: Tue, 14 Oct 2025 14:44:47 +1100 Subject: [PATCH 07/11] :zap: allowing for preprocessing with .gz --- terrier/apps.py | 45 +++++++++++++++++++--------------- terrier/repeatmasker.py | 10 +++++++- tests/test-data/custom.fna.gz | Bin 0 -> 83 bytes 3 files changed, 34 insertions(+), 21 deletions(-) create mode 100644 tests/test-data/custom.fna.gz diff --git a/terrier/apps.py b/terrier/apps.py index 4264506..a61d6f9 100644 --- a/terrier/apps.py +++ b/terrier/apps.py @@ -220,27 +220,32 @@ def checkpoint(self, checkpoint:Path=None) -> str: return checkpoint or "https://github.com/rbturnbull/terrier/releases/download/v0.2.0/terrier-0.2.0.ckpt" def find_fasta_paths(self, files:list[str]) -> list[Path]: - ALLOWABLE_EXTENSIONS = [ - "fasta", - "fa", - "fna", - "fas", - "ref", - ] + base_extensions = {".fa", ".fasta", ".fna", ".fas", ".frn", ".ref"} + + # Function to check if a file matches allowed extensions (including .gz) + def matches_extensions(file: Path): + return ( + file.suffix in base_extensions or + (file.suffix == ".gz" and any(file.stem.endswith(ext) for ext in base_extensions)) + ) + + # Expand the list fasta_paths = [] - for file in files: - file = self.process_location(file) - - if file.is_dir(): - for extension in ALLOWABLE_EXTENSIONS: - fasta_paths += list(file.glob(f'*.{extension}')) - elif file.is_file() and file.exists(): - fasta_paths.append(file) - assert len(fasta_paths), ( - f"No FASTA files found. " - f"Make sure you provide at least one file or a directory of files with an expected FASTA extension:" - f"({', '.join(ALLOWABLE_EXTENSIONS)})" - ) + + # If 'files' is a string or Path, convert it to a list + if isinstance(files, (str,Path)): + fasta_paths = [files] + + for path in files: + path = self.process_location(path) + if path.is_dir(): + # If it's a directory, find all files with the specified extensions + fasta_paths.extend([file for file in path.rglob("*") if matches_extensions(file)]) + else: + # If it's not a directory, add the file to the list + if matches_extensions(path): + fasta_paths.append(path) + return fasta_paths @ta.tool diff --git a/terrier/repeatmasker.py b/terrier/repeatmasker.py index 234a69e..5363b20 100644 --- a/terrier/repeatmasker.py +++ b/terrier/repeatmasker.py @@ -4,6 +4,14 @@ import toml from corgi.seqtree import SeqTree from collections import Counter +import gzip + + +def open_maybe_gz(file:Path): + if file.name.endswith('.gz'): + return gzip.open(file, "rt") + else: + return open(file, "r") def get_verbatim_classification(path:Path, record) -> str: @@ -47,7 +55,7 @@ def create_repeatmasker_seqtree( # Read files count = 0 for file in fasta_paths: - with open(file) as f: + with open_maybe_gz(file) as f: for record in SeqIO.parse(f, "fasta"): partition = count % partitions accession = record.id diff --git a/tests/test-data/custom.fna.gz b/tests/test-data/custom.fna.gz new file mode 100644 index 0000000000000000000000000000000000000000..c3c7f3880c34d3c034d5eee56e735b534f6c55d7 GIT binary patch literal 83 zcmV-Z0IdHXiwFo$#qDSS17me_bZ>1gW^Q2sv-31GGFEo+bJTZCPE1M7<#Kd(58<+l ph)7H+P6J_OpO7GZ=fu*S%-sC4%%aj_E=PCg5CD>$$#6FS007m3B?SNg literal 0 HcmV?d00001 From 7122cbc1420eb953865cbbaf70587105908246f4 Mon Sep 17 00:00:00 2001 From: rturnbull Date: Tue, 14 Oct 2025 15:08:58 +1100 Subject: [PATCH 08/11] :zap: moving find_fasta_paths to corgi --- .gitignore | 3 ++- terrier/apps.py | 29 ----------------------------- 2 files changed, 2 insertions(+), 30 deletions(-) diff --git a/.gitignore b/.gitignore index 9285a09..354c86a 100644 --- a/.gitignore +++ b/.gitignore @@ -185,4 +185,5 @@ oryza-terrier.final.TEs.csv drosophila-terrier.final.TEs.csv lightning_logs/ comparison-test-data/get-rules.ipynb -compare.sh \ No newline at end of file +compare.sh +drosophila.final.TEs.fa.gz \ No newline at end of file diff --git a/terrier/apps.py b/terrier/apps.py index a61d6f9..cb25c5e 100644 --- a/terrier/apps.py +++ b/terrier/apps.py @@ -219,35 +219,6 @@ def get_prediction_probability(row): def checkpoint(self, checkpoint:Path=None) -> str: return checkpoint or "https://github.com/rbturnbull/terrier/releases/download/v0.2.0/terrier-0.2.0.ckpt" - def find_fasta_paths(self, files:list[str]) -> list[Path]: - base_extensions = {".fa", ".fasta", ".fna", ".fas", ".frn", ".ref"} - - # Function to check if a file matches allowed extensions (including .gz) - def matches_extensions(file: Path): - return ( - file.suffix in base_extensions or - (file.suffix == ".gz" and any(file.stem.endswith(ext) for ext in base_extensions)) - ) - - # Expand the list - fasta_paths = [] - - # If 'files' is a string or Path, convert it to a list - if isinstance(files, (str,Path)): - fasta_paths = [files] - - for path in files: - path = self.process_location(path) - if path.is_dir(): - # If it's a directory, find all files with the specified extensions - fasta_paths.extend([file for file in path.rglob("*") if matches_extensions(file)]) - else: - # If it's not a directory, add the file to the list - if matches_extensions(path): - fasta_paths.append(path) - - return fasta_paths - @ta.tool def preprocess( self, From 8310ca3279b336329aad44dd8e021fbc258a53ca Mon Sep 17 00:00:00 2001 From: Robert Turnbull Date: Tue, 14 Oct 2025 15:25:09 +1100 Subject: [PATCH 09/11] :abc: updating docs with more info on training --- README.rst | 12 +++++++++--- docs/reproduction.rst | 16 +++++----------- docs/training.rst | 2 +- 3 files changed, 15 insertions(+), 15 deletions(-) diff --git a/README.rst b/README.rst index 95456ce..5168b90 100644 --- a/README.rst +++ b/README.rst @@ -78,7 +78,7 @@ To run inference on a FASTA file, run this command: .. code-block:: bash - terrier --file INPUT.fa --output-fasta OUTPUT.fa + terrier --input INPUT.fa --output-fasta OUTPUT.fa That will add the classification to after the sequence ID in the `OUTPUT.fa` FASTA file. @@ -86,15 +86,21 @@ If you want to save the probabilities for all classes run this: .. code-block:: bash - terrier --file INPUT.fa --output-csv OUTPUT.csv + terrier --input INPUT.fa --output-csv OUTPUT.csv The columns will be the probability of each classification and the rows correspond to each sequence in ``INPUT.fa``. +You can also use a URL as the input: + +.. code-block:: bash + + terrier --input https://example.com/INPUT.fasta.gz --output-fasta OUTPUT.fa + If you want to output a visualization of the prediction probabilities: .. code-block:: bash - terrier --file INPUT.fa --image-dir OUTPUT-IMAGES/ + terrier --input INPUT.fa --image-dir OUTPUT-IMAGES/ The outputs for the above can be combined together. For more options run diff --git a/docs/reproduction.rst b/docs/reproduction.rst index 2c03567..ab1669c 100644 --- a/docs/reproduction.rst +++ b/docs/reproduction.rst @@ -15,15 +15,12 @@ Fruit Fly Genome Bickmann et al. (2023) provide Transposable Elements (TE) models of a fruit-fly genome. Download it with the following command: -.. code-block:: bash - - wget https://raw.githubusercontent.com/IOB-Muenster/TEclass2/refs/heads/main/tests/Drosophila_melanogaster.fasta Run inference using Terrier like this: .. code-block:: bash - terrier --file Drosophila_melanogaster.fasta \ + terrier --input https://raw.githubusercontent.com/IOB-Muenster/TEclass2/refs/heads/main/tests/drosophila.final.TEs.fa \ --output-csv drosophila-terrier.final.TEs.csv \ --min-length 0 \ --threshold 0 @@ -108,15 +105,12 @@ Rice Genome Bickmann et al. (2023) also provide Transposable Elements (TE) models of a rice genome. Download it with the following command: -.. code-block:: bash - - wget https://raw.githubusercontent.com/IOB-Muenster/TEclass2/refs/heads/main/tests/Oryza_sativa.fasta - Run inference using Terrier like this: .. code-block:: bash - terrier --file Oryza_sativa.fasta --output-csv oryza-terrier.final.TEs.csv --threshold 0 + terrier --intput https://raw.githubusercontent.com/IOB-Muenster/TEclass2/refs/heads/main/tests/oryza.final.TEs.fa \ + --output-csv oryza-terrier.final.TEs.csv --threshold 0 .. note:: @@ -203,7 +197,7 @@ Run inference using Terrier like this: .. code-block:: bash - terrier --file Homo_sapiens.fasta --output-csv Terrier-human.csv --threshold 0 + terrier --input Homo_sapiens.fasta --output-csv Terrier-human.csv --threshold 0 Now evaluate the results with the following command: @@ -243,7 +237,7 @@ Run inference using Terrier like this: .. code-block:: bash - terrier --file Mus_musculus.fasta --output-csv Terrier-mouse.csv --threshold 0 + terrier --input Mus_musculus.fasta --output-csv Terrier-mouse.csv --threshold 0 Now evaluate the results with the following command: diff --git a/docs/training.rst b/docs/training.rst index 59d64d9..759d2cc 100644 --- a/docs/training.rst +++ b/docs/training.rst @@ -7,7 +7,7 @@ After performing the instructions on the :ref:`preprocessing:Preprocessing` page To train Terrier, you will need to use the `terrier-tools` CLI utility. -To use the same hyperparameters as in the main release of Terrier, you can run the following command: +To train with the default settings of Terrier, you can run the following command: .. code-block:: bash From baf2643e06502114b7a0663cb10d0cf337cd77a4 Mon Sep 17 00:00:00 2001 From: Robert Turnbull Date: Tue, 14 Oct 2025 15:25:22 +1100 Subject: [PATCH 10/11] :tv: bumping corgi --- poetry.lock | 8 ++++---- pyproject.toml | 2 +- 2 files changed, 5 insertions(+), 5 deletions(-) diff --git a/poetry.lock b/poetry.lock index 591925a..14e7905 100644 --- a/poetry.lock +++ b/poetry.lock @@ -375,14 +375,14 @@ chardet = ">=3.0.2" [[package]] name = "bio-corgi" -version = "0.5.0a2" +version = "0.5.1" description = "Classifier for ORganelle Genomes Inter alia" optional = false python-versions = "<3.13,>=3.10" groups = ["main"] files = [ - {file = "bio_corgi-0.5.0a2-py3-none-any.whl", hash = "sha256:ead9af52004ba7f581080819456c6aa157b4068683fb7c8551fcad6d81a6e1f1"}, - {file = "bio_corgi-0.5.0a2.tar.gz", hash = "sha256:4c1a6839fef5efff05f475e36b44c34b0939a4c504e5a96f03e1ca2ebdbd6e4a"}, + {file = "bio_corgi-0.5.1-py3-none-any.whl", hash = "sha256:2f3c513e995ffada505c6a63b1db81bad718dc1e4ff39fcb2b9c3a3c95b5edb2"}, + {file = "bio_corgi-0.5.1.tar.gz", hash = "sha256:2d78ab8b6e807c47eb0e63bb8f2a4e065ae9d9a2b3a6a514169391b946e0723b"}, ] [package.dependencies] @@ -5804,4 +5804,4 @@ propcache = ">=0.2.1" [metadata] lock-version = "2.1" python-versions = ">=3.10,<3.13" -content-hash = "30d6d7342864455c1ffe49ced6bbc7b9ce977ff624ad83910d6e771084e45b16" +content-hash = "833ed9506890909ce15aea6a429061fff472694cf28d0fd4f4e0d68755a8ae2c" diff --git a/pyproject.toml b/pyproject.toml index bfe9f74..d286251 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -26,7 +26,7 @@ python = ">=3.10,<3.13" numpy = "<2.0.0" pandas = "<=2.2.0" h5py = ">=3.8.0" -bio-corgi = ">=0.5.0a2" +bio-corgi = ">=0.5.1" pyfastx = ">=1.1.0" toml = ">=0.10.2" From 0e31dadc3e4fefe4a5e44326b02cbde3a4e30be8 Mon Sep 17 00:00:00 2001 From: Robert Turnbull Date: Tue, 14 Oct 2025 15:26:28 +1100 Subject: [PATCH 11/11] :tv: bumping version' --- pyproject.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pyproject.toml b/pyproject.toml index d286251..b533720 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -1,6 +1,6 @@ [tool.poetry] name = "bio-terrier" -version = "0.3.3" +version = "0.3.4" description = "Transposable Element Repeat Result classifIER" authors = ["Robert Turnbull "] license = "Apache-2.0"