diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml new file mode 100644 index 0000000..dceaed8 --- /dev/null +++ b/.github/workflows/ci.yml @@ -0,0 +1,38 @@ +name: Smoke Test + +on: + push: + branches: [master, add_tests] + pull_request: + branches: [master] + +jobs: + smoke-test: + runs-on: ubuntu-latest + steps: + - uses: actions/checkout@v4 + + - name: Set up Python + uses: actions/setup-python@v5 + with: + python-version: "3.11" + + - name: Install dependencies + run: | + python -m pip install --upgrade pip + pip install biopython numpy gffutils pyfaidx + + - name: Verify Python scripts parse without errors + run: | + find pipelines -name '*.py' | while read f; do + echo "Checking $f ..." + python -c "import py_compile; py_compile.compile('$f', doraise=True)" + done + + - name: Verify Nextflow config is valid + uses: nf-core/setup-nextflow@v2 + with: + version: "latest-stable" + + - name: Nextflow smoke test + run: nextflow run main.nf -help diff --git a/.github/workflows/tests.yml b/.github/workflows/tests.yml new file mode 100644 index 0000000..f662ce6 --- /dev/null +++ b/.github/workflows/tests.yml @@ -0,0 +1,26 @@ +name: Unit Tests + +on: + push: + branches: [master, add_tests] + pull_request: + branches: [master] + +jobs: + unit-tests: + runs-on: ubuntu-latest + steps: + - uses: actions/checkout@v4 + + - name: Set up Python + uses: actions/setup-python@v5 + with: + python-version: "3.11" + + - name: Install dependencies + run: | + python -m pip install --upgrade pip + pip install biopython numpy gffutils pyfaidx pytest pytest-cov + + - name: Run unit tests with coverage + run: pytest tests/ -v --cov=pipelines --cov-report=term-missing --cov-report=xml diff --git a/README.md b/README.md index 5f2f146..165f564 100644 --- a/README.md +++ b/README.md @@ -1,22 +1,117 @@ +[![Smoke Test](https://github.com/CAMI-challenge/CAMISIM/actions/workflows/ci.yml/badge.svg)](https://github.com/CAMI-challenge/CAMISIM/actions/workflows/ci.yml) +[![Unit Tests](https://github.com/CAMI-challenge/CAMISIM/actions/workflows/tests.yml/badge.svg)](https://github.com/CAMI-challenge/CAMISIM/actions/workflows/tests.yml) + # CAMISIM -CAMISIM is a software to model abundance distributions of microbial communities and to simulate corresponding shotgun metagenome datasets.\ -It was mainly developed for the [Critical Assessment of Metagenome Annotation (CAMI)](http://microbiome-cosi.org/cami) challenge, but should be suitable for general use. Please don't hesitate to [open a new issue](https://github.com/CAMI-challenge/CAMISIM/issues) if you run into problems or need help. +CAMISIM simulates realistic shotgun metagenome and metatranscriptome datasets from microbial communities. It models abundance distributions, generates reads using established simulators, and produces gold standard assemblies and binning references for benchmarking. + +Originally developed for the [Critical Assessment of Metagenome Interpretation (CAMI)](http://microbiome-cosi.org/cami) challenge, CAMISIM is suitable for general use. Please [open an issue](https://github.com/CAMI-challenge/CAMISIM/issues) if you run into problems. + +## What's new in v2.0 + +CAMISIM 2.0 is a complete rewrite using [Nextflow](https://www.nextflow.io/) (DSL2), replacing the previous Python-based pipeline. Key changes: + +- Nextflow-managed workflow with support for local, SLURM, and cloud executors +- Conda/Mamba-based dependency management +- Metatranscriptomic simulation pipeline (new) +- Configurable read simulators: ART (Illumina), NanoSim3 (Nanopore), WGSIM, pbsim3 (PacBio) + +The previous standalone Python version is available via the git tag `1.31-final` but will not receive updates. Use `convert_config.py` to migrate v1 config files to v2 format. + +## Quick start + +### Requirements + +- [Nextflow](https://www.nextflow.io/) (DSL2) +- [Conda](https://docs.conda.io/) or [Mamba](https://mamba.readthedocs.io/) (for automatic environment management) + +### Run a metagenomic simulation + +```bash +nextflow run main.nf --pipeline metagenomic +``` + +### Run a metatranscriptomic simulation + +```bash +nextflow run main.nf --pipeline metatranscriptomic +``` + +### Use a custom config + +```bash +nextflow run main.nf --pipeline metagenomic --config path/to/your.config +``` + +Both pipelines ship with default example data under `nextflow_defaults/` so you can do a test run out of the box. + +## Key parameters + +| Parameter | Default | Description | +|-----------|---------|-------------| +| `size` | `0.05` | Sample size in Gbp | +| `type` | `art` | Read simulator (`art`, `nanosim3`, `wgsim`, `pbsim3`) | +| `number_of_samples` | `2` | Number of samples to simulate | +| `gsa` | `true` | Generate gold standard assembly | +| `anonymization` | `true` | Anonymize read and contig IDs | +| `seed` | `632741178` | Random seed for reproducibility | +| `biom_profile` | `""` | BIOM profile for community design (optional) | + +See [Configuration File Options](https://github.com/CAMI-challenge/CAMISIM/wiki/Configuration-File-Options) for the full parameter reference. + +## Input files + +| File | Format | Description | +|------|--------|-------------| +| `genome_locations.tsv` | TSV | `genome_IDpath/to/genome.fa` | +| `meta_data.tsv` | TSV with header | `genome_IDOTUNCBI_IDnovelty_category` | +| `distribution_N.txt` | TSV | `genome_IDabundance_fraction` (one per sample) | +| `gene_annotation_locations.tsv` | TSV | `genome_IDpath/to/annotations.gff3` (metatranscriptomic only) | + +## Project structure + +``` +CAMISIM/ + main.nf # Entry point + nextflow.config # Global config (selects pipeline) + pipelines/ + metagenomic/ # Metagenomic simulation pipeline + metagenomic.nf # Main workflow + sample_wise_simulation.nf # Per-sample read simulation + read_simulators/ # ART, NanoSim3, WGSIM modules + scripts/ # Python utilities + config/ # Default parameters + metatranscriptomic/ # Metatranscriptomic simulation pipeline + metatranscriptomic.nf # Main workflow (+ gene expression) + read_simulators/ # ART, NanoSim3, pbsim3 modules + scripts/ # Python utilities + config/ # Default parameters + shared/ # Shared modules (both pipelines) + anonymization.nf # Read/contig anonymization + binning.nf # Read-to-genome mapping + distribution.nf # Abundance normalization + scripts/ # Shared Python utilities + nextflow_defaults/ # Example input data for test runs + tools/ # Bundled simulators and reference data + tests/ # Unit tests (pytest) +``` + +## Documentation + +- [User manual](https://github.com/CAMI-challenge/CAMISIM/wiki/User-manual) +- [Configuration File Options](https://github.com/CAMI-challenge/CAMISIM/wiki/Configuration-File-Options) +- [File Formats](https://github.com/CAMI-challenge/CAMISIM/wiki/File-Formats) + +## Citation + +If you use CAMISIM, please cite: -### CAMISIM 2.0 -CAMISIM received a major update to version 2.0, using nextflow. You can still use the old python standalone using the tag 1.31-final, it will not receive updates however.\ -You can use the script "convert_config.py" to convert your CAMISIM1 config file to CAMISIM2 (no guarantees for correctness, please check yourself).\ -This version has been tested, but if you encounter any unforeseen difficulties or differences from what you expect your simulation to look like, please raise an Issue. +> Fritz\*, Hofmann\*, et al. (2019). **CAMISIM: Simulating metagenomes and microbial communities.** *Microbiome*, 7:17. doi:[10.1186/s40168-019-0633-6](https://doi.org/10.1186/s40168-019-0633-6) -### Documentation -* [User manual](https://github.com/CAMI-challenge/CAMISIM/wiki/User-manual) -* [Configuration File Options](https://github.com/CAMI-challenge/CAMISIM/wiki/Configuration-File-Options) -* [File Formats](https://github.com/CAMI-challenge/CAMISIM/wiki/File-Formats) +You may also cite the CAMI benchmark paper: -### Citation +> Sczyrba\*, Hofmann\*, Belmann\*, et al. (2017). **Critical Assessment of Metagenome Interpretation — a benchmark of metagenomics software.** *Nature Methods*, 14(11):1063–1071. doi:[10.1038/nmeth.4458](https://doi.org/10.1038/nmeth.4458) -If you use CAMISIM, please cite the publication at *Microbiome*: -* Fritz*, Hofmann*, et al. (2019). **CAMISIM: Simulating metagenomes and microbial communities.** *Microbiome*, 2019, 7:17. doi:[10.1186/s40168-019-0633-6](https://doi.org/10.1186/s40168-019-0633-6) +## License -A part of CAMISIM's functionality was also described in the CAMI manuscript, thus you may also cite: -* Sczyrba*, Hofmann*, Belmann*, et al. (2017). **Critical Assessment of Metagenome Interpretation—a benchmark of metagenomics software.** *Nature Methods*, 14, 11:1063–1071. doi:[10.1038/nmeth.4458](https://doi.org/10.1038/nmeth.4458) +Apache License 2.0. See [LICENSE.txt](LICENSE.txt). diff --git a/tests/conftest.py b/tests/conftest.py new file mode 100644 index 0000000..0aeea27 --- /dev/null +++ b/tests/conftest.py @@ -0,0 +1,10 @@ +"""Shared fixtures and path setup for CAMISIM tests.""" + +import sys +import os + +# Add script directories to Python path so we can import modules +ROOT = os.path.dirname(os.path.dirname(os.path.abspath(__file__))) +sys.path.insert(0, os.path.join(ROOT, "pipelines", "shared", "scripts")) +sys.path.insert(0, os.path.join(ROOT, "pipelines", "metagenomic", "scripts")) +sys.path.insert(0, os.path.join(ROOT, "pipelines", "metatranscriptomic", "scripts")) diff --git a/tests/test_anonymizer.py b/tests/test_anonymizer.py new file mode 100644 index 0000000..cb96d68 --- /dev/null +++ b/tests/test_anonymizer.py @@ -0,0 +1,120 @@ +"""Tests for pipelines/shared/scripts/anonymizer.py — sequence anonymization.""" + +import io +import pytest +from anonymizer import Anonymizer + + +@pytest.fixture +def anon(): + return Anonymizer() + + +def _make_fasta(*seqs): + """Create a FASTA-formatted StringIO from (id, seq) tuples.""" + lines = [] + for sid, seq in seqs: + lines.append(f">{sid}\n{seq}\n") + return io.StringIO("".join(lines)) + + +class TestAnonymizeSequences: + def test_basic_anonymization(self, anon): + inp = _make_fasta(("seq1", "ACGT"), ("seq2", "TTTT")) + out = io.StringIO() + mapping = io.StringIO() + anon.anonymize_sequences(mapping, input_stream=inp, output_stream=out, file_format="fasta") + + mapping.seek(0) + lines = mapping.read().strip().split("\n") + assert len(lines) == 2 + assert lines[0].startswith("seq1\t") + assert lines[1].startswith("seq2\t") + + def test_prefix(self, anon): + inp = _make_fasta(("seq1", "ACGT")) + out = io.StringIO() + mapping = io.StringIO() + anon.anonymize_sequences(mapping, input_stream=inp, output_stream=out, sequence_prefix="S0R", file_format="fasta") + + mapping.seek(0) + assert "S0R0" in mapping.read() + + def test_sequential_ids(self, anon): + inp = _make_fasta(("a", "A"), ("b", "C"), ("c", "G")) + out = io.StringIO() + mapping = io.StringIO() + anon.anonymize_sequences(mapping, input_stream=inp, output_stream=out, file_format="fasta") + + mapping.seek(0) + ids = [line.split("\t")[1] for line in mapping.read().strip().split("\n")] + assert ids == ["0", "1", "2"] + + def test_output_contains_sequences(self, anon): + inp = _make_fasta(("seq1", "ACGT")) + out = io.StringIO() + mapping = io.StringIO() + anon.anonymize_sequences(mapping, input_stream=inp, output_stream=out, file_format="fasta") + + out.seek(0) + content = out.read() + assert "ACGT" in content + + def test_invalid_format_raises(self, anon): + with pytest.raises(AssertionError): + anon.anonymize_sequences(io.StringIO(), file_format="bam") + + +class TestAnonymizeSequencePairs: + def test_paired_ids(self, anon): + inp = _make_fasta(("r1/1", "ACGT"), ("r1/2", "TTTT"), ("r2/1", "GGGG"), ("r2/2", "CCCC")) + out = io.StringIO() + mapping = io.StringIO() + anon.anonymize_sequence_pairs(mapping, input_stream=inp, output_stream=out, file_format="fasta") + + mapping.seek(0) + new_ids = [line.split("\t")[1] for line in mapping.read().strip().split("\n")] + # pair 0: 0/1, 0/2; pair 1: 1/1, 1/2 + assert new_ids == ["0/1", "0/2", "1/1", "1/2"] + + def test_prefix_with_pairs(self, anon): + inp = _make_fasta(("a", "A"), ("b", "C")) + out = io.StringIO() + mapping = io.StringIO() + anon.anonymize_sequence_pairs(mapping, input_stream=inp, output_stream=out, sequence_prefix="X", file_format="fasta") + + mapping.seek(0) + new_ids = [line.split("\t")[1] for line in mapping.read().strip().split("\n")] + assert new_ids == ["X0/1", "X0/2"] + + def test_six_sequences_three_pairs(self, anon): + inp = _make_fasta( + ("r1f", "AAAA"), ("r1r", "TTTT"), + ("r2f", "GGGG"), ("r2r", "CCCC"), + ("r3f", "ACGT"), ("r3r", "TGCA"), + ) + out = io.StringIO() + mapping = io.StringIO() + anon.anonymize_sequence_pairs(mapping, input_stream=inp, output_stream=out, file_format="fasta") + + mapping.seek(0) + new_ids = [line.split("\t")[1] for line in mapping.read().strip().split("\n")] + assert new_ids == ["0/1", "0/2", "1/1", "1/2", "2/1", "2/2"] + + +class TestAnonymizerEdgeCases: + def test_empty_input(self, anon): + inp = io.StringIO("") + out = io.StringIO() + mapping = io.StringIO() + anon.anonymize_sequences(mapping, input_stream=inp, output_stream=out, file_format="fasta") + assert out.getvalue() == "" + assert mapping.getvalue() == "" + + def test_original_ids_in_mapping(self, anon): + inp = _make_fasta(("my_original_id", "ACGT")) + out = io.StringIO() + mapping = io.StringIO() + anon.anonymize_sequences(mapping, input_stream=inp, output_stream=out, file_format="fasta") + mapping.seek(0) + assert "my_original_id" in mapping.read() diff --git a/tests/test_build_ncbi_taxonomy.py b/tests/test_build_ncbi_taxonomy.py new file mode 100644 index 0000000..040553c --- /dev/null +++ b/tests/test_build_ncbi_taxonomy.py @@ -0,0 +1,336 @@ +"""Tests for pipelines/metagenomic/scripts/build_ncbi_taxonomy.py — taxonomy utilities.""" + +import io +import sys +import os +import pytest + +sys.path.insert(0, os.path.join(os.path.dirname(__file__), "..", "pipelines", "metagenomic", "scripts")) +import build_ncbi_taxonomy as tax + + +@pytest.fixture(autouse=True) +def _reset_globals(): + """Reset global state before each test.""" + tax.TAXID_TO_NAME.clear() + tax.TAXID_TO_PARENT_TAXID.clear() + tax.TAXID_TO_RANK.clear() + tax.TAXID_OLD_TO_TAXID_NEW.clear() + tax.IS_LEGACY = False + yield + tax.TAXID_TO_NAME.clear() + tax.TAXID_TO_PARENT_TAXID.clear() + tax.TAXID_TO_RANK.clear() + tax.TAXID_OLD_TO_TAXID_NEW.clear() + tax.IS_LEGACY = False + + +def _setup_simple_taxonomy(): + """Set up a minimal taxonomy: 1 -> 2 (superkingdom) -> 100 (genus) -> 200 (species).""" + tax.TAXID_TO_PARENT_TAXID.update({ + "1": "1", + "2": "1", + "100": "2", + "200": "100", + }) + tax.TAXID_TO_RANK.update({ + "1": "no rank", + "2": "superkingdom", + "100": "genus", + "200": "species", + }) + tax.TAXID_TO_NAME.update({ + "1": "root", + "2": "Bacteria", + "100": "Escherichia", + "200": "Escherichia coli", + }) + + +class TestGetUpdatedTaxid: + def test_known_taxid(self): + tax.TAXID_TO_RANK["12345"] = "species" + assert tax.get_updated_taxid("12345") == "12345" + + def test_merged_taxid(self): + tax.TAXID_TO_RANK["999"] = "species" + tax.TAXID_OLD_TO_TAXID_NEW["111"] = "999" + assert tax.get_updated_taxid("111") == "999" + + def test_special_taxid_other_entries(self): + assert tax.get_updated_taxid(tax.OTHER_ENTRIES_TAXID) == tax.OTHER_ENTRIES_TAXID + + def test_special_taxid_plasmid(self): + assert tax.get_updated_taxid(tax.PLASMID_SPECIES_TAXID) == tax.PLASMID_SPECIES_TAXID + + def test_unknown_taxid_raises(self): + with pytest.raises(ValueError): + tax.get_updated_taxid("9999999") + + +class TestGetScientificName: + def test_known_name(self): + tax.TAXID_TO_RANK["200"] = "species" + tax.TAXID_TO_NAME["200"] = "Escherichia coli" + assert tax.get_scientific_name("200") == "Escherichia coli" + + def test_merged_taxid_lookup(self): + tax.TAXID_TO_RANK["200"] = "species" + tax.TAXID_TO_NAME["200"] = "Escherichia coli" + tax.TAXID_OLD_TO_TAXID_NEW["111"] = "200" + assert tax.get_scientific_name("111") == "Escherichia coli" + + def test_unknown_raises(self): + with pytest.raises(ValueError): + tax.get_scientific_name("9999999") + + +class TestIsDescendantOf: + def test_direct_descendant(self): + _setup_simple_taxonomy() + assert tax.is_descendant_of("200", "100") is True + + def test_transitive_descendant(self): + _setup_simple_taxonomy() + assert tax.is_descendant_of("200", "2") is True + + def test_not_descendant(self): + _setup_simple_taxonomy() + assert tax.is_descendant_of("2", "200") is False + + def test_self_is_descendant(self): + _setup_simple_taxonomy() + assert tax.is_descendant_of("200", "200") is True + + def test_unknown_taxid(self): + assert tax.is_descendant_of("999", "1") is False + + +class TestGetOtherEntriesLineage: + def test_basic_lineage(self): + lineage = tax.get_other_entries_lineage("45202") + assert lineage[0] == tax.OTHER_ENTRIES_TAXID + assert lineage[8] == "45202" + assert lineage[9] == "" + assert len(lineage) == len(tax.OTHER_RANKS) + + def test_with_strain(self): + lineage = tax.get_other_entries_lineage("45202", taxid_if_strain="99999") + assert lineage[9] == "99999" + + def test_empty_middle_ranks(self): + lineage = tax.get_other_entries_lineage("45202") + for i in range(1, 8): + assert lineage[i] == "" + + +class TestGetLineageOfLegalRanks: + def test_legacy_ranks(self): + tax.IS_LEGACY = True + _setup_simple_taxonomy() + lineage = tax.get_lineage_of_legal_ranks("200") + assert len(lineage) == len(tax.LEGACY_RANKS) + # species should be filled + species_idx = tax.LEGACY_RANKS.index("species") + assert lineage[species_idx] == "200" + # genus should be filled + genus_idx = tax.LEGACY_RANKS.index("genus") + assert lineage[genus_idx] == "100" + + def test_missing_ranks_get_default(self): + tax.IS_LEGACY = True + _setup_simple_taxonomy() + lineage = tax.get_lineage_of_legal_ranks("200", default_value="") + # phylum etc. not in our mini taxonomy → default + phylum_idx = tax.LEGACY_RANKS.index("phylum") + assert lineage[phylum_idx] == "" + + +class TestParseStream: + def test_basic_tsv(self): + lines = ["a\tb\tc\n", "d\te\tf\n"] + rows = list(tax.parse_stream(lines)) + assert rows == [["a", "b", "c"], ["d", "e", "f"]] + + def test_skips_comments(self): + lines = ["# comment\n", "a\tb\n"] + rows = list(tax.parse_stream(lines)) + assert len(rows) == 1 + + def test_skips_empty_lines(self): + lines = ["\n", "a\tb\n", "\n"] + rows = list(tax.parse_stream(lines)) + assert len(rows) == 1 + + def test_strips_newlines(self): + lines = ["a\tb\n"] + rows = list(tax.parse_stream(lines)) + assert rows[0] == ["a", "b"] + + +class TestStreamTpHeader: + def test_header_format(self): + out = io.StringIO() + tax.stream_tp_header(out, "sample_0") + content = out.getvalue() + assert "@SampleID:sample_0" in content + assert f"@Version:{tax.TAXONOMIC_PROFILE_VERSION}" in content + + +class TestGetOtherEntriesLineageNames: + def test_plasmid_species(self): + tax.TAXID_TO_NAME[tax.PLASMID_SPECIES_TAXID] = tax.PLASMID_SPECIES_NAME + names = tax.get_other_entries_lineage_names(tax.PLASMID_SPECIES_TAXID) + assert names[0] == tax.OTHER_ENTRIES_NAME + assert names[8] == tax.PLASMID_SPECIES_NAME + + def test_with_strain(self): + tax.TAXID_TO_NAME["12345"] = "Some species" + names = tax.get_other_entries_lineage_names("12345", taxid_if_strain="99999") + assert "strain" in names[9] + + +class TestReadNamesFile: + def test_parses_scientific_names(self, tmp_path): + names_file = tmp_path / "names.dmp" + names_file.write_text( + "65\t|\tHerpetosiphon aurantiacus\t|\t\t|\tscientific name\t|\n" + "65\t|\tH. aurantiacus\t|\t\t|\tsynonym\t|\n" + "100\t|\tEscherichia\t|\t\t|\tscientific name\t|\n" + ) + tax.read_names_file(str(names_file)) + assert tax.TAXID_TO_NAME["65"] == "Herpetosiphon aurantiacus" + assert tax.TAXID_TO_NAME["100"] == "Escherichia" + # synonym should not overwrite + assert tax.TAXID_TO_NAME["65"] != "H. aurantiacus" + + def test_injects_special_nodes(self, tmp_path): + names_file = tmp_path / "names.dmp" + names_file.write_text("") + tax.read_names_file(str(names_file)) + assert tax.TAXID_TO_NAME[tax.OTHER_ENTRIES_TAXID] == tax.OTHER_ENTRIES_NAME + assert tax.TAXID_TO_NAME[tax.PLASMID_SPECIES_TAXID] == tax.PLASMID_SPECIES_NAME + + +class TestBuildNcbiTaxonomy: + def test_parses_nodes(self, tmp_path): + nodes_file = tmp_path / "nodes.dmp" + nodes_file.write_text( + "1\t|\t1\t|\tno rank\t|\n" + "2\t|\t1\t|\tsuperkingdom\t|\n" + "100\t|\t2\t|\tgenus\t|\n" + ) + tax.build_ncbi_taxonomy(str(nodes_file)) + assert tax.TAXID_TO_PARENT_TAXID["2"] == "1" + assert tax.TAXID_TO_RANK["100"] == "genus" + assert tax.IS_LEGACY is True # no "acellular root" + + def test_detects_new_taxonomy(self, tmp_path): + nodes_file = tmp_path / "nodes.dmp" + nodes_file.write_text( + "1\t|\t1\t|\tno rank\t|\n" + "999\t|\t1\t|\tacellular root\t|\n" + ) + tax.build_ncbi_taxonomy(str(nodes_file)) + assert tax.IS_LEGACY is False + + +class TestReadMergedFile: + def test_parses_merges(self, tmp_path): + merged_file = tmp_path / "merged.dmp" + merged_file.write_text("111\t|\t999\t|\n222\t|\t888\t|\n") + tax.read_merged_file(str(merged_file)) + assert tax.TAXID_OLD_TO_TAXID_NEW["111"] == "999" + assert tax.TAXID_OLD_TO_TAXID_NEW["222"] == "888" + + +class TestParseFile: + def test_reads_tsv(self, tmp_path): + f = tmp_path / "data.tsv" + f.write_text("a\tb\n# comment\nc\td\n") + rows = list(tax.parse_file(str(f))) + assert rows == [["a", "b"], ["c", "d"]] + + def test_empty_file(self, tmp_path): + f = tmp_path / "empty.tsv" + f.write_text("") + rows = list(tax.parse_file(str(f))) + assert rows == [] + + +class TestGetPercentByRankByTaxid: + def test_legacy_mode(self): + tax.IS_LEGACY = True + _setup_simple_taxonomy() + # lineage aligned to LEGACY_RANKS: [superkingdom, phylum, class, order, family, genus, species, strain] + lineage = [None] * len(tax.LEGACY_RANKS) + lineage[tax.LEGACY_RANKS.index("species")] = "200" + lineage[tax.LEGACY_RANKS.index("genus")] = "100" + lineage[tax.LEGACY_RANKS.index("superkingdom")] = "2" + + result = tax.get_percent_by_rank_by_taxid( + genome_id_to_lineage={"g1": lineage}, + genome_id_to_percent={"g1": 0.75}, + genome_id_is_plasmid={"g1": False}, + ) + assert "legacy" in result + assert result["legacy"]["species"]["200"] == 0.75 + assert result["legacy"]["genus"]["100"] == 0.75 + + def test_two_genomes_same_species(self): + tax.IS_LEGACY = True + _setup_simple_taxonomy() + lineage = [None] * len(tax.LEGACY_RANKS) + lineage[tax.LEGACY_RANKS.index("species")] = "200" + lineage[tax.LEGACY_RANKS.index("superkingdom")] = "2" + + result = tax.get_percent_by_rank_by_taxid( + genome_id_to_lineage={"g1": lineage, "g2": lineage}, + genome_id_to_percent={"g1": 0.3, "g2": 0.7}, + genome_id_is_plasmid={"g1": False, "g2": False}, + ) + assert abs(result["legacy"]["species"]["200"] - 1.0) < 1e-9 + + +class TestGetGenomeIdToLineage: + def test_regular_genome(self): + tax.IS_LEGACY = True + _setup_simple_taxonomy() + result = tax.get_genome_id_to_lineage( + list_of_genome_id=["g1"], + genome_id_to_taxid={"g1": "200"}, + strain_id_to_genome_id={}, + genome_id_to_strain_id={}, + genome_id_is_plasmid={"g1": False}, + ) + assert "g1" in result + lineage = result["g1"] + species_idx = tax.LEGACY_RANKS.index("species") + assert lineage[species_idx] == "200" + # strain should be auto-assigned + strain_idx = tax.LEGACY_RANKS.index("strain") + assert lineage[strain_idx] is not None + + def test_plasmid_genome(self): + _setup_simple_taxonomy() + result = tax.get_genome_id_to_lineage( + list_of_genome_id=["p1"], + genome_id_to_taxid={"p1": "200"}, + strain_id_to_genome_id={}, + genome_id_to_strain_id={}, + genome_id_is_plasmid={"p1": True}, + ) + lineage = result["p1"] + assert lineage[0] == tax.OTHER_ENTRIES_TAXID + assert lineage[8] == tax.PLASMID_SPECIES_TAXID + + def test_empty_taxid_raises(self): + _setup_simple_taxonomy() + with pytest.raises(KeyError, match="has no taxid"): + tax.get_genome_id_to_lineage( + list_of_genome_id=["g1"], + genome_id_to_taxid={"g1": ""}, + strain_id_to_genome_id={}, + genome_id_to_strain_id={}, + genome_id_is_plasmid={"g1": False}, + ) diff --git a/tests/test_clean_up_sequences.py b/tests/test_clean_up_sequences.py new file mode 100644 index 0000000..7e561a5 --- /dev/null +++ b/tests/test_clean_up_sequences.py @@ -0,0 +1,28 @@ +"""Tests for pipelines/metagenomic/scripts/clean_up_sequences.py""" + +import pytest +from clean_up_sequences import get_new_name + + +class TestGetNewName: + def test_no_collision(self): + assert get_new_name("seq1", set()) == "seq1" + + def test_single_collision(self): + names = {"seq1"} + result = get_new_name("seq1", names) + assert result == "seq1_0" + + def test_multiple_collisions(self): + names = {"seq1", "seq1_0", "seq1_1"} + result = get_new_name("seq1", names) + assert result == "seq1_2" + + def test_no_collision_with_other_names(self): + names = {"seq2", "seq3"} + assert get_new_name("seq1", names) == "seq1" + + def test_result_not_in_existing_set(self): + names = {"contig", "contig_0"} + result = get_new_name("contig", names) + assert result not in names diff --git a/tests/test_community_distribution.py b/tests/test_community_distribution.py new file mode 100644 index 0000000..4e68bc9 --- /dev/null +++ b/tests/test_community_distribution.py @@ -0,0 +1,479 @@ +"""Tests for pipelines/shared/scripts/get_community_distribution.py""" + +import io +import random +import pytest +import get_community_distribution as cd + + +# --------------------------------------------------------------------------- +# lt_zero +# --------------------------------------------------------------------------- +class TestLtZero: + def test_positive_value_unchanged(self): + assert cd.lt_zero(5.0) == 5.0 + + def test_zero_returns_floor(self): + assert cd.lt_zero(0) == 0.001 + + def test_negative_returns_floor(self): + assert cd.lt_zero(-3.5) == 0.001 + + def test_small_positive_unchanged(self): + assert cd.lt_zero(0.0001) == 0.0001 + + +# --------------------------------------------------------------------------- +# Boolean helpers +# --------------------------------------------------------------------------- +class TestIsBooleanState: + @pytest.mark.parametrize("word", ["yes", "true", "on", "y", "t"]) + def test_truthy_words(self, word): + assert cd.is_boolean_state(word) is True + + @pytest.mark.parametrize("word", ["no", "false", "off", "n", "f"]) + def test_falsy_words(self, word): + assert cd.is_boolean_state(word) is True + + @pytest.mark.parametrize("word", ["maybe", "1", "0", "YES", "True", ""]) + def test_non_boolean_words(self, word): + assert cd.is_boolean_state(word) is False + + +class TestGetBooleanState: + @pytest.mark.parametrize("word,expected", [ + ("yes", True), ("true", True), ("on", True), ("y", True), ("t", True), + ("no", False), ("false", False), ("off", False), ("n", False), ("f", False), + ]) + def test_valid_words(self, word, expected): + assert cd.get_boolean_state(word) is expected + + def test_invalid_word_raises(self): + with pytest.raises(AssertionError): + cd.get_boolean_state("maybe") + + +class TestStrToBool: + def test_true(self): + assert cd.str_to_bool("True") is True + + def test_false(self): + assert cd.str_to_bool("False") is False + + def test_invalid_raises(self): + with pytest.raises(ValueError): + cd.str_to_bool("true") + + def test_empty_raises(self): + with pytest.raises(ValueError): + cd.str_to_bool("") + + +# --------------------------------------------------------------------------- +# validate_number +# --------------------------------------------------------------------------- +class TestValidateNumber: + def test_valid_number_no_bounds(self): + assert cd.validate_number(5) is True + + def test_below_minimum(self): + assert cd.validate_number(1, minimum=5) is False + + def test_above_maximum(self): + assert cd.validate_number(10, maximum=5) is False + + def test_within_range(self): + assert cd.validate_number(5, minimum=1, maximum=10) is True + + def test_zero_excluded(self): + assert cd.validate_number(0, zero=False) is False + + def test_zero_allowed_by_default(self): + assert cd.validate_number(0) is True + + def test_float_valid(self): + assert cd.validate_number(3.14, minimum=1.0, maximum=4.0) is True + + def test_non_number_raises(self): + with pytest.raises(AssertionError): + cd.validate_number("five") + + def test_at_exact_minimum(self): + assert cd.validate_number(5, minimum=5) is True + + def test_at_exact_maximum(self): + assert cd.validate_number(5, maximum=5) is True + + def test_negative_number(self): + assert cd.validate_number(-3, minimum=-10, maximum=0) is True + + +# --------------------------------------------------------------------------- +# get_initial_list +# --------------------------------------------------------------------------- +class TestGetInitialList: + def test_dimensions(self): + result = cd.get_initial_list(3, 2) + assert len(result) == 3 + assert all(len(row) == 2 for row in result) + + def test_all_zeros(self): + result = cd.get_initial_list(2, 3) + for row in result: + assert all(v == 0.0 for v in row) + + def test_independent_rows(self): + result = cd.get_initial_list(2, 2) + result[0][0] = 1.0 + assert result[1][0] == 0.0 + + def test_non_int_raises(self): + with pytest.raises(AssertionError): + cd.get_initial_list(2.0, 3) + + def test_single_element(self): + result = cd.get_initial_list(1, 1) + assert result == [[0.0]] + + +# --------------------------------------------------------------------------- +# has_unique_columns +# --------------------------------------------------------------------------- +class TestHasUniqueColumns: + def test_unique(self): + assert cd.has_unique_columns(["a", "b", "c"]) is True + + def test_duplicates(self): + assert cd.has_unique_columns(["a", "b", "a"]) is False + + def test_empty(self): + assert cd.has_unique_columns([]) is True + + def test_single(self): + assert cd.has_unique_columns(["a"]) is True + + +# --------------------------------------------------------------------------- +# get_distribution_file_paths +# --------------------------------------------------------------------------- +class TestGetDistributionFilePaths: + def test_correct_count(self): + result = cd.get_distribution_file_paths("/out", 3) + assert len(result) == 3 + + def test_correct_format(self): + result = cd.get_distribution_file_paths("/out", 2) + assert result[0].endswith("distribution_0.txt") + assert result[1].endswith("distribution_1.txt") + + def test_directory_prefix(self): + result = cd.get_distribution_file_paths("/my/dir", 1) + assert result[0].startswith("/my/dir") + + def test_zero_samples(self): + result = cd.get_distribution_file_paths("/out", 0) + assert result == [] + + +# --------------------------------------------------------------------------- +# get_compression_type +# --------------------------------------------------------------------------- +class TestGetCompressionType: + def test_gz(self): + assert cd.get_compression_type("data.fastq.gz") == "gz" + + def test_bz2(self): + assert cd.get_compression_type("data.tar.bz2") == "bz2" + + def test_no_compression(self): + assert cd.get_compression_type("data.txt") is None + + def test_unknown_extension(self): + assert cd.get_compression_type("data.xyz") is None + + +# --------------------------------------------------------------------------- +# has_column / get_map +# --------------------------------------------------------------------------- +class TestHasColumn: + def test_existing_column(self): + table = {"name": ["a", "b"], "id": [1, 2]} + assert cd.has_column(table, "name") is True + + def test_missing_column(self): + table = {"name": ["a", "b"]} + assert cd.has_column(table, "id") is False + + def test_int_key(self): + table = {0: ["a"], 1: ["b"]} + assert cd.has_column(table, 0) is True + + +class TestGetMap: + def test_basic_mapping(self): + table = {"key": ["a", "b", "c"], "val": [1, 2, 3]} + result = cd.get_map(table, "key", "val") + assert result == {"a": 1, "b": 2, "c": 3} + + def test_duplicate_key_raises(self): + table = {"key": ["a", "a"], "val": [1, 2]} + with pytest.raises(KeyError): + cd.get_map(table, "key", "val", unique_key=True) + + def test_missing_column_raises(self): + table = {"key": ["a"]} + with pytest.raises(AssertionError): + cd.get_map(table, "key", "val") + + def test_empty_table(self): + table = {"key": []} + result = cd.get_map(table, "key", "key") + assert result == {} + + def test_same_column_for_key_and_value(self): + # get_map requires len(metadata_table) >= 2 columns to produce results + table = {"col": ["x", "y", "z"], "other": ["a", "b", "c"]} + result = cd.get_map(table, "col", "col") + assert result == {"x": "x", "y": "y", "z": "z"} + + +# --------------------------------------------------------------------------- +# Distribution functions (seeded for reproducibility) +# --------------------------------------------------------------------------- +class TestAddInitialLogDistribution: + def test_fills_first_sample(self): + pop = cd.get_initial_list(5, 3) + random.seed(42) + cd.add_initial_log_distribution(pop, 1.0, 1.0) + for row in pop: + assert row[0] > 0 # lognormal is always positive + assert row[1] == 0.0 # other samples untouched + assert row[2] == 0.0 + + def test_reproducible_with_seed(self): + pop1 = cd.get_initial_list(3, 1) + random.seed(99) + cd.add_initial_log_distribution(pop1, 2.0, 0.5) + + pop2 = cd.get_initial_list(3, 1) + random.seed(99) + cd.add_initial_log_distribution(pop2, 2.0, 0.5) + assert pop1 == pop2 + + +class TestAddReplicates: + def test_fills_remaining_samples(self): + pop = cd.get_initial_list(3, 4) + random.seed(42) + cd.add_initial_log_distribution(pop, 1.0, 1.0) + cd.add_replicates(pop, 0, 0.1) + for row in pop: + assert all(v > 0 for v in row) + + def test_first_sample_unchanged(self): + pop = cd.get_initial_list(2, 3) + random.seed(42) + cd.add_initial_log_distribution(pop, 1.0, 1.0) + first_vals = [row[0] for row in pop] + cd.add_replicates(pop, 0, 0.1) + for i, row in enumerate(pop): + assert row[0] == first_vals[i] + + +class TestAddTimeseriesGauss: + def test_extinction_propagates(self): + pop = [[0.0, 0.0, 0.0]] + cd.add_timeseries_gauss(pop, 0, 1.0) + assert pop[0][1] == 0.0 + assert pop[0][2] == 0.0 + + def test_positive_values_produce_positive(self): + pop = [[5.0, 0.0, 0.0]] + random.seed(42) + cd.add_timeseries_gauss(pop, 0, 0.001) + assert all(v > 0 for v in pop[0]) + + +class TestAddTimeseriesLognorm: + def test_fills_sequentially(self): + pop = cd.get_initial_list(2, 3) + random.seed(42) + cd.add_initial_log_distribution(pop, 1.0, 1.0) + cd.add_timeseries_lognorm(pop, 1.0, 1.0) + for row in pop: + assert all(v > 0 for v in row) + + +class TestAddDifferential: + def test_fills_independently(self): + pop = cd.get_initial_list(2, 3) + random.seed(42) + cd.add_initial_log_distribution(pop, 1.0, 1.0) + cd.add_differential(pop, 1.0, 1.0) + for row in pop: + assert all(v > 0 for v in row[1:]) + + +class TestRandomDistributionToRelativeAbundance: + def test_sums_to_one(self): + pop = [[3.0, 6.0], [7.0, 4.0]] + cd.random_distribution_to_relative_abundance(pop) + for sample_idx in range(2): + total = sum(row[sample_idx] for row in pop) + assert abs(total - 1.0) < 1e-9 + + def test_proportions_correct(self): + pop = [[2.0], [8.0]] + cd.random_distribution_to_relative_abundance(pop) + assert abs(pop[0][0] - 0.2) < 1e-9 + assert abs(pop[1][0] - 0.8) < 1e-9 + + def test_single_genome(self): + pop = [[5.0, 3.0]] + cd.random_distribution_to_relative_abundance(pop) + assert pop[0][0] == 1.0 + assert pop[0][1] == 1.0 + + +# --------------------------------------------------------------------------- +# get_lists_of_distributions (integration of the above) +# --------------------------------------------------------------------------- +class TestGetListsOfDistributions: + @pytest.mark.parametrize("modus", [ + "replicates", "timeseries_normal", "timeseries_lognormal", "differential", + ]) + def test_all_modes_produce_valid_distributions(self, modus): + random.seed(42) + result = cd.get_lists_of_distributions( + size_of_population=5, + number_of_samples=3, + modus=modus, + log_mu=1.0, + log_sigma=1.0, + gauss_mu=0.0, + gauss_sigma=0.5, + ) + assert len(result) == 5 + assert all(len(row) == 3 for row in result) + # each sample column should sum to ~1.0 + for s in range(3): + total = sum(row[s] for row in result) + assert abs(total - 1.0) < 1e-6 + + def test_deterministic_with_seed(self): + random.seed(123) + r1 = cd.get_lists_of_distributions(3, 2, "replicates", 1.0, 1.0, 0.0, 0.5) + random.seed(123) + r2 = cd.get_lists_of_distributions(3, 2, "replicates", 1.0, 1.0, 0.0, 0.5) + assert r1 == r2 + + +# --------------------------------------------------------------------------- +# write_distribution_file +# --------------------------------------------------------------------------- +class TestWriteDistributionFile: + def test_basic_output(self): + out = io.StringIO() + abundances = {"genome_A": [0.3, 0.7], "genome_B": [0.6, 0.4]} + cd.write_distribution_file(out, abundances, 1) + content = out.getvalue() + assert "genome_A\t0.3\n" in content + assert "genome_B\t0.6\n" in content + + def test_second_sample(self): + out = io.StringIO() + abundances = {"g1": [0.1, 0.9]} + cd.write_distribution_file(out, abundances, 2) + assert "g1\t0.9\n" in out.getvalue() + + def test_all_genomes_written(self): + out = io.StringIO() + abundances = {f"g{i}": [float(i)] for i in range(10)} + cd.write_distribution_file(out, abundances, 1) + lines = [l for l in out.getvalue().strip().split("\n") if l] + assert len(lines) == 10 + + +# --------------------------------------------------------------------------- +# read (file-based TSV reader) +# --------------------------------------------------------------------------- +class TestRead: + def test_basic_tsv(self, tmp_path): + f = tmp_path / "data.tsv" + f.write_text("a\t1\nb\t2\nc\t3\n") + result = cd.read(str(f)) + assert result[0] == ["a", "b", "c"] + assert result[1] == ["1", "2", "3"] + + def test_skips_comments(self, tmp_path): + f = tmp_path / "data.tsv" + f.write_text("# header\na\t1\nb\t2\n") + result = cd.read(str(f)) + assert result[0] == ["a", "b"] + + def test_skips_empty_lines(self, tmp_path): + f = tmp_path / "data.tsv" + f.write_text("a\t1\n\nb\t2\n") + result = cd.read(str(f)) + assert result[0] == ["a", "b"] + + def test_custom_separator(self, tmp_path): + f = tmp_path / "data.csv" + f.write_text("a,1\nb,2\n") + result = cd.read(str(f), separator=",") + assert result[0] == ["a", "b"] + assert result[1] == ["1", "2"] + + def test_custom_comment_char(self, tmp_path): + f = tmp_path / "data.tsv" + f.write_text("% comment\na\t1\n") + result = cd.read(str(f), comment_line="%") + assert result[0] == ["a"] + + +# --------------------------------------------------------------------------- +# openy (file opener with compression support) +# --------------------------------------------------------------------------- +class TestOpeny: + def test_open_plain_text(self, tmp_path): + f = tmp_path / "test.txt" + f.write_text("hello") + fh = cd.openy(str(f)) + assert fh.read() == "hello" + fh.close() + + def test_invalid_mode_raises(self): + with pytest.raises(AssertionError): + cd.openy("test.txt", mode="x") + + def test_write_and_read_gz(self, tmp_path): + f = str(tmp_path / "test.gz") + fh = cd.openy(f, mode="w", compression_type="gz") + fh.write(b"compressed data") + fh.close() + fh = cd.openy(f, mode="r") + assert b"compressed data" in fh.read() + fh.close() + + +# --------------------------------------------------------------------------- +# get_genome_id_to_path_map +# --------------------------------------------------------------------------- +class TestGetGenomeIdToPathMap: + def test_basic_mapping(self, tmp_path): + f = tmp_path / "genome_locations.tsv" + f.write_text("genome1\t/path/to/g1.fa\ngenome2\t/path/to/g2.fa\n") + result = cd.get_genome_id_to_path_map(str(f)) + assert result == {"genome1": "/path/to/g1.fa", "genome2": "/path/to/g2.fa"} + + def test_empty_file(self, tmp_path): + f = tmp_path / "empty.tsv" + f.write_text("") + result = cd.get_genome_id_to_path_map(str(f)) + assert result == {} + + def test_single_genome(self, tmp_path): + f = tmp_path / "single.tsv" + f.write_text("gA\t/genomes/a.fasta\n") + result = cd.get_genome_id_to_path_map(str(f)) + assert result == {"gA": "/genomes/a.fasta"} diff --git a/tests/test_convert_config.py b/tests/test_convert_config.py new file mode 100644 index 0000000..8931842 --- /dev/null +++ b/tests/test_convert_config.py @@ -0,0 +1,121 @@ +"""Tests for convert_config.py — v1 INI to v2 Nextflow config conversion.""" + +import sys +import os +import pytest + +sys.path.insert(0, os.path.join(os.path.dirname(__file__), "..")) +from convert_config import convert_to_nextflow_config + + +class TestConvertToNextflowConfig: + def test_output_directory_renamed(self, tmp_path): + ini = tmp_path / "test.ini" + ini.write_text("[Main]\noutput_directory = /my/output\n") + out = tmp_path / "out.config" + convert_to_nextflow_config(str(ini), str(out)) + content = out.read_text() + assert 'outdir = "/my/output"' in content + + def test_gsa_true(self, tmp_path): + ini = tmp_path / "test.ini" + ini.write_text("[Main]\ngsa = True\n") + out = tmp_path / "out.config" + convert_to_nextflow_config(str(ini), str(out)) + assert "gsa = true" in out.read_text() + + def test_gsa_false(self, tmp_path): + ini = tmp_path / "test.ini" + ini.write_text("[Main]\ngsa = False\n") + out = tmp_path / "out.config" + convert_to_nextflow_config(str(ini), str(out)) + assert "gsa = false" in out.read_text() + + def test_anonymous_renamed_to_anonymization(self, tmp_path): + ini = tmp_path / "test.ini" + ini.write_text("[Main]\nanonymous = True\n") + out = tmp_path / "out.config" + convert_to_nextflow_config(str(ini), str(out)) + assert "anonymization = true" in out.read_text() + + def test_type_quoted(self, tmp_path): + ini = tmp_path / "test.ini" + ini.write_text("[Main]\ntype = nanosim3\n") + out = tmp_path / "out.config" + convert_to_nextflow_config(str(ini), str(out)) + assert 'type = "nanosim3"' in out.read_text() + + def test_skipped_keys(self, tmp_path): + ini = tmp_path / "test.ini" + ini.write_text("[Main]\nphase = 0\nmax_processors = 8\ncompress = True\nsamtools = /usr/bin/samtools\n") + out = tmp_path / "out.config" + convert_to_nextflow_config(str(ini), str(out)) + content = out.read_text() + assert "phase" not in content.split("biom_profile")[0] # before defaults + assert "max_processors" not in content.split("biom_profile")[0] + assert "compress" not in content.split("biom_profile")[0] + + def test_metadata_renamed(self, tmp_path): + ini = tmp_path / "test.ini" + ini.write_text("[Main]\nmetadata = /path/to/meta.tsv\n") + out = tmp_path / "out.config" + convert_to_nextflow_config(str(ini), str(out)) + assert 'metadata_file = "/path/to/meta.tsv"' in out.read_text() + + def test_id_to_genome_file_renamed(self, tmp_path): + ini = tmp_path / "test.ini" + ini.write_text("[Main]\nid_to_genome_file = /path/to/genomes.tsv\n") + out = tmp_path / "out.config" + convert_to_nextflow_config(str(ini), str(out)) + assert 'genome_locations_file = "/path/to/genomes.tsv"' in out.read_text() + + def test_error_profiles_renamed(self, tmp_path): + ini = tmp_path / "test.ini" + ini.write_text("[Main]\nerror_profiles = /path/to/profiles\n") + out = tmp_path / "out.config" + convert_to_nextflow_config(str(ini), str(out)) + assert 'base_profile_name = "/path/to/profiles"' in out.read_text() + + def test_fragment_size_sd_renamed(self, tmp_path): + ini = tmp_path / "test.ini" + ini.write_text("[Main]\nfragment_size_standard_deviation = 27\n") + out = tmp_path / "out.config" + convert_to_nextflow_config(str(ini), str(out)) + assert "fragment_size_sd = 27" in out.read_text() + + def test_pooled_gsa_true(self, tmp_path): + ini = tmp_path / "test.ini" + ini.write_text("[Main]\npooled_gsa = True\n") + out = tmp_path / "out.config" + convert_to_nextflow_config(str(ini), str(out)) + assert "pooled_gsa = true" in out.read_text() + + def test_pooled_gsa_false(self, tmp_path): + ini = tmp_path / "test.ini" + ini.write_text("[Main]\npooled_gsa = False\n") + out = tmp_path / "out.config" + convert_to_nextflow_config(str(ini), str(out)) + assert 'pooled_gsa = "[]"' in out.read_text() + + def test_view_renamed_to_verbose(self, tmp_path): + ini = tmp_path / "test.ini" + ini.write_text("[Main]\nview = True\n") + out = tmp_path / "out.config" + convert_to_nextflow_config(str(ini), str(out)) + assert "verbose = true" in out.read_text() + + def test_defaults_appended(self, tmp_path): + ini = tmp_path / "test.ini" + ini.write_text("[Main]\nsize = 0.1\n") + out = tmp_path / "out.config" + convert_to_nextflow_config(str(ini), str(out)) + content = out.read_text() + assert 'biom_profile=""' in content + assert "conda.enabled = true" in content + + def test_passthrough_key(self, tmp_path): + ini = tmp_path / "test.ini" + ini.write_text("[Main]\nnumber_of_samples = 5\n") + out = tmp_path / "out.config" + convert_to_nextflow_config(str(ini), str(out)) + assert "number_of_samples = 5" in out.read_text() diff --git a/tests/test_gene_abundance.py b/tests/test_gene_abundance.py new file mode 100644 index 0000000..74b8246 --- /dev/null +++ b/tests/test_gene_abundance.py @@ -0,0 +1,94 @@ +"""Tests for pipelines/metatranscriptomic/scripts/get_gene_abundance.py.""" + +import random +import pytest + +import sys, os +sys.path.insert(0, os.path.join(os.path.dirname(__file__), "..", "pipelines", "metatranscriptomic", "scripts")) +from get_gene_abundance import GeneAbundace + + +@pytest.fixture +def ga(): + return GeneAbundace() + + +class TestLtZero: + def test_positive(self, ga): + assert ga.lt_zero(5.0) == 5.0 + + def test_zero(self, ga): + assert ga.lt_zero(0) == 0.001 + + def test_negative(self, ga): + assert ga.lt_zero(-1.0) == 0.001 + + +class TestAddInitialLogDistribution: + def test_fills_specified_sample(self, ga): + abundances = {"g1": [0.0, 0.0], "g2": [0.0, 0.0]} + random.seed(42) + ga.add_initial_log_distribution(abundances, 1.0, 1.0, 0) + for val in abundances.values(): + assert val[0] > 0 + assert val[1] == 0.0 + + def test_fills_second_sample(self, ga): + abundances = {"g1": [0.0, 0.0]} + random.seed(42) + ga.add_initial_log_distribution(abundances, 1.0, 1.0, 1) + assert abundances["g1"][0] == 0.0 + assert abundances["g1"][1] > 0 + + def test_reproducible(self, ga): + a1 = {"g1": [0.0], "g2": [0.0]} + random.seed(99) + ga.add_initial_log_distribution(a1, 2.0, 0.5, 0) + + a2 = {"g1": [0.0], "g2": [0.0]} + random.seed(99) + ga.add_initial_log_distribution(a2, 2.0, 0.5, 0) + assert a1 == a2 + + +class TestAddTimeseriesGauss: + def test_extinction_propagates(self, ga): + abundances = {"g1": [0.0, 0.0, 0.0]} + ga.add_timeseries_gauss(abundances, 0, 1.0) + assert abundances["g1"][1] == 0.0 + assert abundances["g1"][2] == 0.0 + + def test_positive_chain(self, ga): + abundances = {"g1": [5.0, 0.0, 0.0]} + random.seed(42) + ga.add_timeseries_gauss(abundances, 0, 0.001) + assert all(v > 0 for v in abundances["g1"]) + + +class TestAddGenewiseLogDistribution: + def test_fills_from_initial(self, ga): + abundances = {"g1": [2.0, 0.0, 0.0]} + random.seed(42) + ga.add_genewise_log_distribution(abundances, 0.5) + assert all(v > 0 for v in abundances["g1"]) + + +class TestRandomDistributionToRelativeAbundance: + def test_sums_to_one(self, ga): + abundances = {"g1": [3.0, 6.0], "g2": [7.0, 4.0]} + ga.random_distribution_to_relative_abundance(abundances, 2) + for s in range(2): + total = sum(v[s] for v in abundances.values()) + assert abs(total - 1.0) < 1e-9 + + def test_proportions(self, ga): + abundances = {"g1": [2.0], "g2": [8.0]} + ga.random_distribution_to_relative_abundance(abundances, 1) + assert abs(abundances["g1"][0] - 0.2) < 1e-9 + assert abs(abundances["g2"][0] - 0.8) < 1e-9 + + def test_single_gene(self, ga): + abundances = {"g1": [5.0, 3.0]} + ga.random_distribution_to_relative_abundance(abundances, 2) + assert abundances["g1"][0] == 1.0 + assert abundances["g1"][1] == 1.0 diff --git a/tests/test_get_genomes.py b/tests/test_get_genomes.py new file mode 100644 index 0000000..69d8852 --- /dev/null +++ b/tests/test_get_genomes.py @@ -0,0 +1,93 @@ +"""Tests for pipelines/metagenomic/scripts/get_genomes.py — pure utility functions.""" + +import sys +import os +import pytest + +sys.path.insert(0, os.path.join(os.path.dirname(__file__), "..", "pipelines", "metagenomic", "scripts")) + +pytest.importorskip("biom", reason="biom not installed") +from get_genomes import sort_by_abundance, read_genomes_list, str2bool, truncated_geometric + + +class TestSortByAbundance: + def test_sorts_descending(self): + profile = { + "otu1": (["k__Bacteria"], [1.0, 2.0]), + "otu2": (["k__Archaea"], [10.0, 20.0]), + "otu3": (["k__Fungi"], [5.0, 5.0]), + } + result = sort_by_abundance(profile) + assert result == ["otu2", "otu3", "otu1"] + + def test_single_otu(self): + profile = {"otu1": ([], [3.0])} + assert sort_by_abundance(profile) == ["otu1"] + + def test_equal_abundance(self): + profile = { + "a": ([], [5.0]), + "b": ([], [5.0]), + } + result = sort_by_abundance(profile) + assert len(result) == 2 + assert set(result) == {"a", "b"} + + +class TestReadGenomesList: + def test_basic_parsing(self, tmp_path): + f = tmp_path / "genomes.tsv" + f.write_text("12345\tE. coli\tftp://example.com/path\n") + result, total = read_genomes_list(str(f)) + assert total == 1 + assert "12345" in result + name, paths, novelty = result["12345"] + assert name == "E. coli" + assert "http://example.com/path" in paths[0] # ftp replaced with http + assert novelty == "known_strain" + + def test_multiple_genomes_same_taxid(self, tmp_path): + f = tmp_path / "genomes.tsv" + f.write_text("100\tSpeciesA\tftp://a/path1\n100\tSpeciesA\tftp://a/path2\n") + result, total = read_genomes_list(str(f)) + assert total == 2 + assert len(result["100"][1]) == 2 + + def test_with_additional_file(self, tmp_path): + main = tmp_path / "genomes.tsv" + main.write_text("100\tSpeciesA\tftp://a/path\n") + add = tmp_path / "additional.tsv" + add.write_text("200\tSpeciesB\t/local/genome.fa\tnovel_strain\n") + result, total = read_genomes_list(str(main), str(add)) + assert total == 2 + assert "200" in result + assert result["200"][2] == "novel_strain" + + +class TestStr2Bool: + @pytest.mark.parametrize("val", ["1", "true", "True", "TRUE", "yes", "y", "Y"]) + def test_truthy(self, val): + assert str2bool(val) is True + + @pytest.mark.parametrize("val", ["0", "false", "False", "no", "n", ""]) + def test_falsy(self, val): + assert str2bool(val) is False + + +class TestTruncatedGeometric: + def test_within_bounds(self): + for _ in range(50): + val = truncated_geometric(0.5, 1, 5) + assert 1 <= val <= 5 + + def test_single_value_range(self): + val = truncated_geometric(0.5, 3, 3) + assert val == 3 + + def test_invalid_params_raises(self): + with pytest.raises(ValueError): + truncated_geometric(0, 1, 5) + + def test_invalid_range_raises(self): + with pytest.raises(ValueError): + truncated_geometric(0.5, 10, 1) diff --git a/tests/test_goldstandardfileformat.py b/tests/test_goldstandardfileformat.py new file mode 100644 index 0000000..a9ba0dc --- /dev/null +++ b/tests/test_goldstandardfileformat.py @@ -0,0 +1,222 @@ +"""Tests for pipelines/shared/scripts/goldstandardfileformat.py — TSV parsing.""" + +import io +import pytest +from goldstandardfileformat import GoldStandardFileFormat + + +@pytest.fixture +def gsf(): + return GoldStandardFileFormat() + + +class TestParseColumnNames: + def test_basic_header(self, gsf): + stream = io.StringIO("name\tid\tvalue\nrow1\n") + cols = gsf._parse_column_names(stream, "\t") + assert cols == ["name", "id", "value"] + + def test_duplicate_columns_raises(self, gsf): + stream = io.StringIO("name\tname\n") + with pytest.raises(AssertionError): + gsf._parse_column_names(stream, "\t") + + +class TestProcessLines: + def test_without_column_names(self, gsf): + stream = io.StringIO("a\tb\nc\td\n") + result = gsf.process_lines(stream, {}, ["#"], "\t", column_names=False) + assert result[0] == ["a", "c"] + assert result[1] == ["b", "d"] + + def test_with_column_names(self, gsf): + stream = io.StringIO("name\tid\nalice\t1\nbob\t2\n") + result = gsf.process_lines(stream, {}, ["#"], "\t", column_names=True) + assert result["name"] == ["alice", "bob"] + assert result["id"] == ["1", "2"] + + def test_skips_comments(self, gsf): + stream = io.StringIO("# comment\na\tb\n") + result = gsf.process_lines(stream, {}, ["#"], "\t", column_names=False) + assert result[0] == ["a"] + + def test_skips_empty_lines(self, gsf): + stream = io.StringIO("a\tb\n\nc\td\n") + result = gsf.process_lines(stream, {}, ["#"], "\t", column_names=False) + assert result[0] == ["a", "c"] + + def test_bad_column_count_raises(self, gsf): + stream = io.StringIO("name\tid\nalice\t1\nbob\n") + with pytest.raises(ValueError, match="Bad number of values"): + gsf.process_lines(stream, {}, ["#"], "\t", column_names=True) + + +class TestRead: + def test_read_from_stream(self, gsf): + stream = io.TextIOWrapper(io.BytesIO(b"x\ty\n1\t2\n3\t4\n")) + result = gsf.read(stream, separator="\t", column_names=True) + assert result["x"] == ["1", "3"] + assert result["y"] == ["2", "4"] + + def test_read_without_headers(self, gsf, tmp_path): + f = tmp_path / "data.tsv" + f.write_text("a\tb\nc\td\n") + result = gsf.read(str(f), separator="\t") + assert result[0] == ["a", "c"] + + +class TestHasColumn: + def test_existing(self, gsf): + assert gsf.has_column({"a": [1]}, "a") is True + + def test_missing(self, gsf): + assert gsf.has_column({"a": [1]}, "b") is False + + def test_int_key(self, gsf): + assert gsf.has_column({0: ["x"]}, 0) is True + + +class TestGetMap: + def test_basic(self, gsf): + table = {"k": ["a", "b"], "v": [1, 2]} + assert gsf.get_map(table, "k", "v") == {"a": 1, "b": 2} + + def test_duplicate_key_raises(self, gsf): + table = {"k": ["a", "a"], "v": [1, 2]} + with pytest.raises(KeyError): + gsf.get_map(table, "k", "v", unique_key=True) + + def test_missing_column_raises(self, gsf): + with pytest.raises(AssertionError): + gsf.get_map({"k": []}, "k", "missing") + + +class TestWriteGsReadMapping: + def test_basic_mapping(self, gsf): + out = io.StringIO() + anon_to_read = {"A0": "seq1.1-0", "A1": "seq1.1-1"} + seq_to_genome = {"seq1.1": "genome1"} + genome_to_tax = {"genome1": "12345"} + gsf.write_gs_read_mapping( + out, anon_to_read, seq_to_genome, genome_to_tax, + nanosim_real_fastq=False, wgsim=False, + ) + content = out.getvalue() + assert "#anonymous_read_id" in content + assert "A0\tgenome1\t12345\tseq1.1-0" in content + + def test_wgsim_mode(self, gsf): + out = io.StringIO() + anon_to_read = {"A0": "seq1_100_200_0:0:0_0:0:0_1/1"} + seq_to_genome = {"seq1": "g1"} + genome_to_tax = {"g1": "999"} + gsf.write_gs_read_mapping( + out, anon_to_read, seq_to_genome, genome_to_tax, + nanosim_real_fastq=False, wgsim=True, + ) + content = out.getvalue() + assert "g1\t999" in content + + def test_missing_dash_raises(self, gsf): + out = io.StringIO() + anon_to_read = {"A0": "readwithoutdash"} + seq_to_genome = {} + genome_to_tax = {} + with pytest.raises(ValueError, match="missing '-'"): + gsf.write_gs_read_mapping( + out, anon_to_read, seq_to_genome, genome_to_tax, + nanosim_real_fastq=False, wgsim=False, + ) + + def test_missing_underscore_wgsim_raises(self, gsf): + out = io.StringIO() + anon_to_read = {"A0": "readwithoutunderscore"} + seq_to_genome = {} + genome_to_tax = {} + with pytest.raises(ValueError, match="missing '_'"): + gsf.write_gs_read_mapping( + out, anon_to_read, seq_to_genome, genome_to_tax, + nanosim_real_fastq=False, wgsim=True, + ) + + +class TestWriteGsaContigMapping: + def test_basic_contig_mapping(self, gsf): + out = io.StringIO() + seq_to_anon = {"seq1_from_100_to_200_stuff": "contig_0"} + seq_positions = {"seq1": [50, 100, 150, 200, 250]} + seq_to_genome = {"seq1": "g1"} + genome_to_tax = {"g1": "555"} + gsf.write_gsa_contig_mapping( + out, seq_to_anon, seq_positions, seq_to_genome, genome_to_tax, + ) + content = out.getvalue() + assert "#anonymous_contig_id" in content + assert "contig_0\tg1\t555\tseq1" in content + + def test_bad_sequence_id_raises(self, gsf): + out = io.StringIO() + seq_to_anon = {"unknown_from_0_to_100_x": "c0"} + seq_positions = {"unknown": [50]} + seq_to_genome = {} # missing + genome_to_tax = {} + with pytest.raises(KeyError): + gsf.write_gsa_contig_mapping( + out, seq_to_anon, seq_positions, seq_to_genome, genome_to_tax, + ) + + +class TestGetDictUniqueIdToGenomeFilePath: + def test_basic_mapping(self, gsf, tmp_path): + f = tmp_path / "mapping.tsv" + f.write_text("genome1\t/path/to/g1.fa\ngenome2\t/path/to/g2.fa\n") + result = gsf.get_dict_unique_id_to_genome_file_path(str(f)) + assert result == {"genome1": "/path/to/g1.fa", "genome2": "/path/to/g2.fa"} + + +class TestGetDictGenomeIdToTaxId: + def test_basic_metadata(self, gsf, tmp_path): + f = tmp_path / "metadata.tsv" + f.write_text("genome_ID\tOTU\tNCBI_ID\tnovelty\ng1\totu1\t12345\tknown\ng2\totu2\t67890\tnovel\n") + result = gsf.get_dict_genome_id_to_tax_id(str(f)) + assert result == {"g1": "12345", "g2": "67890"} + + +class TestGetDictAnonymousToOriginalId: + def test_reverse_mapping(self, gsf, tmp_path): + f = tmp_path / "anon_map.tsv" + f.write_text("original_seq1\tanon_0\noriginal_seq2\tanon_1\n") + result = gsf.get_dict_anonymous_to_original_id(str(f)) + assert result == {"anon_0": "original_seq1", "anon_1": "original_seq2"} + + +class TestGetDictSequenceNameToAnonymous: + def test_forward_mapping(self, gsf, tmp_path): + f = tmp_path / "seq_map.tsv" + f.write_text("original_seq1\tanon_0\noriginal_seq2\tanon_1\n") + result = gsf.get_dict_sequence_name_to_anonymous(str(f)) + assert result == {"original_seq1": "anon_0", "original_seq2": "anon_1"} + + +class TestGetDictSequenceNameToPositions: + def test_basic_positions(self, gsf, tmp_path): + f = tmp_path / "positions.tsv" + f.write_text("seq1-0\t100\nseq1-1\t200\nseq2-0\t50\n") + result = gsf.get_dict_sequence_name_to_positions([str(f)]) + assert "seq1" in result + assert sorted(result["seq1"]) == [100, 200] + assert result["seq2"] == [50] + + def test_wgsim_format(self, gsf, tmp_path): + f = tmp_path / "positions.tsv" + f.write_text("seq1_100_200_0:0:0_0:0:0_1\t150\n") + result = gsf.get_dict_sequence_name_to_positions([str(f)], wgsim=True) + assert "seq1" in result + + def test_multiple_files(self, gsf, tmp_path): + f1 = tmp_path / "pos1.tsv" + f1.write_text("seq1-0\t100\n") + f2 = tmp_path / "pos2.tsv" + f2.write_text("seq1-1\t200\n") + result = gsf.get_dict_sequence_name_to_positions([str(f1), str(f2)]) + assert sorted(result["seq1"]) == [100, 200] diff --git a/tests/test_maf_converter.py b/tests/test_maf_converter.py new file mode 100644 index 0000000..ba57971 --- /dev/null +++ b/tests/test_maf_converter.py @@ -0,0 +1,36 @@ +"""Tests for pipelines/metatranscriptomic/scripts/maf_converter.py — CIGAR creation.""" + +import sys +import os +import pytest + +sys.path.insert(0, os.path.join(os.path.dirname(__file__), "..", "pipelines", "metatranscriptomic", "scripts")) +from maf_converter import cigar_code_creation + + +class TestCigarCodeCreation: + def test_all_matches(self): + assert cigar_code_creation("ACGT", "ACGT") == "4M" + + def test_mismatch_counted_as_match(self): + assert cigar_code_creation("ACGT", "AGGT") == "4M" + + def test_insertion_gap_in_ref(self): + # gap in reference = insertion in read + assert cigar_code_creation("-ACG", "TACG") == "1I3M" + + def test_deletion_gap_in_read(self): + # gap in read = deletion from reference + assert cigar_code_creation("ACGT", "A-GT") == "1M1D2M" + + def test_multiple_insertions(self): + assert cigar_code_creation("--ACGT", "TTACGT") == "2I4M" + + def test_multiple_deletions(self): + assert cigar_code_creation("ACGT", "A--T") == "1M2D1M" + + def test_single_base_match(self): + assert cigar_code_creation("A", "A") == "1M" + + def test_empty_strings(self): + assert cigar_code_creation("", "") == "" diff --git a/tests/test_pick_random_strains.py b/tests/test_pick_random_strains.py new file mode 100644 index 0000000..0c4fb84 --- /dev/null +++ b/tests/test_pick_random_strains.py @@ -0,0 +1,34 @@ +"""Tests for pipelines/metagenomic/scripts/pick_random_strains.py — helper utilities.""" + +import sys +import os +import pytest + +sys.path.insert(0, os.path.join(os.path.dirname(__file__), "..", "pipelines", "metagenomic", "scripts")) +from pick_random_strains import get_empty_row + + +class TestGetEmptyRow: + def test_returns_dict_by_default(self): + result = get_empty_row(["a", "b", "c"]) + assert result == {"a": "", "b": "", "c": ""} + + def test_custom_default_value(self): + result = get_empty_row(["x", "y"], default_value="NA") + assert result == {"x": "NA", "y": "NA"} + + def test_as_list(self): + result = get_empty_row(["a", "b", "c"], as_list=True) + assert result == ["", "", ""] + + def test_as_list_with_default(self): + result = get_empty_row(["a", "b"], default_value="0", as_list=True) + assert result == ["0", "0"] + + def test_empty_columns(self): + assert get_empty_row([]) == {} + assert get_empty_row([], as_list=True) == [] + + def test_non_string_default_raises(self): + with pytest.raises(AssertionError): + get_empty_row(["a"], default_value=0) diff --git a/tests/test_prepare_strain_simulation.py b/tests/test_prepare_strain_simulation.py new file mode 100644 index 0000000..15923ef --- /dev/null +++ b/tests/test_prepare_strain_simulation.py @@ -0,0 +1,74 @@ +"""Tests for pipelines/metagenomic/scripts/prepare_strain_simulation.py — strain distribution.""" + +import random +import sys +import os +import pytest + +sys.path.insert(0, os.path.join(os.path.dirname(__file__), "..", "pipelines", "metagenomic", "scripts")) +from prepare_strain_simulation import PrepareStrainSimulation + + +@pytest.fixture +def pss(): + return PrepareStrainSimulation() + + +class TestInit: + def test_init_without_seed(self, pss): + pss.init(seed=None) # should not raise + + def test_init_with_seed(self, pss): + pss.init(seed=42) # should not raise + + def test_seed_reproducibility(self, pss): + pss.init(seed=123) + val1 = random.random() + pss.init(seed=123) + val2 = random.random() + assert val1 == val2 + + +class TestGetGenomeAmountsGeometric: + def test_sum_equals_max(self, pss): + pss.init(seed=42) + amounts = pss.get_genome_amounts_geometric(0.5, 20) + assert sum(amounts) == 20 + + def test_all_positive(self, pss): + pss.init(seed=42) + amounts = pss.get_genome_amounts_geometric(0.5, 10) + assert all(a > 0 for a in amounts) + + def test_high_probability_many_singletons(self, pss): + pss.init(seed=42) + amounts = pss.get_genome_amounts_geometric(0.99, 100) + singletons = sum(1 for a in amounts if a == 1) + assert singletons > len(amounts) * 0.5 + + def test_small_max(self, pss): + pss.init(seed=42) + amounts = pss.get_genome_amounts_geometric(0.5, 1) + assert sum(amounts) == 1 + + +class TestGetGenomeAmountsGeometricFix: + def test_sum_equals_max(self, pss): + pss.init(seed=42) + amounts = pss.get_genome_amounts_geometric_fix(5, 20) + assert sum(amounts) == 20 + + def test_correct_number_of_entries(self, pss): + pss.init(seed=42) + amounts = pss.get_genome_amounts_geometric_fix(5, 20) + assert len(amounts) == 5 + + def test_all_at_least_one(self, pss): + pss.init(seed=42) + amounts = pss.get_genome_amounts_geometric_fix(3, 10) + assert all(a >= 1 for a in amounts) + + def test_equal_when_max_equals_real(self, pss): + pss.init(seed=42) + amounts = pss.get_genome_amounts_geometric_fix(5, 5) + assert amounts == [1, 1, 1, 1, 1] diff --git a/tests/test_sam_from_reads.py b/tests/test_sam_from_reads.py new file mode 100644 index 0000000..6c025e6 --- /dev/null +++ b/tests/test_sam_from_reads.py @@ -0,0 +1,130 @@ +"""Tests for pipelines/shared/scripts/sam_from_reads.py.""" + +import io +import os +import sys +import tempfile + +import pytest + +sys.path.insert(0, os.path.join(os.path.dirname(__file__), "..", "pipelines", "shared", "scripts")) +from sam_from_reads import SamFromReads + + +@pytest.fixture +def sam(): + """Create a SamFromReads instance for testing.""" + return SamFromReads() + + +# --------------------------------------------------------------------------- +# get_cigar_length +# --------------------------------------------------------------------------- +class TestGetCigarLength: + def test_simple_match(self, sam): + assert sam.get_cigar_length("10M") == 10 + + def test_multi_digit(self, sam): + assert sam.get_cigar_length("150M") == 150 + + def test_match_and_insertion(self, sam): + assert sam.get_cigar_length("5M3I2M") == 10 + + def test_deletion_excluded(self, sam): + assert sam.get_cigar_length("5M2D3M") == 8 + + def test_complex_cigar(self, sam): + assert sam.get_cigar_length("10M2D5I3M") == 18 + + def test_single_base(self, sam): + assert sam.get_cigar_length("1M") == 1 + + def test_soft_clip(self, sam): + assert sam.get_cigar_length("5S10M5S") == 20 + + def test_only_insertion(self, sam): + assert sam.get_cigar_length("5I") == 5 + + def test_only_deletion(self, sam): + assert sam.get_cigar_length("5D") == 0 + + def test_large_numbers(self, sam): + assert sam.get_cigar_length("1000M50I20D30M") == 1080 + + +# --------------------------------------------------------------------------- +# get_cigars_nanosim +# --------------------------------------------------------------------------- +class TestGetCigarsNanosim: + def _write_error_profile(self, tmp_path, lines): + ep = tmp_path / "error_profile.txt" + ep.write_text("\n".join(lines) + "\n") + return str(ep) + + def test_insertion_only(self, sam, tmp_path): + lines = [ + "Seq\tPos\tType\tLen\tRef\tQuery", + "refA_100_aligned_0_F_0_500_0\t10\tins\t3\tACG\tACGTTT", + ] + ep = self._write_error_profile(tmp_path, lines) + cigars = sam.get_cigars_nanosim(ep) + assert "refA-100-0" in cigars + cigar_str, _ = cigars["refA-100-0"] + assert "10M" in cigar_str + assert "3I" in cigar_str + + def test_deletion_only(self, sam, tmp_path): + lines = [ + "Seq\tPos\tType\tLen\tRef\tQuery", + "refA_100_aligned_0_F_0_500_0\t10\tdel\t2\tAC\t--", + ] + ep = self._write_error_profile(tmp_path, lines) + cigars = sam.get_cigars_nanosim(ep) + cigar_str, _ = cigars["refA-100-0"] + assert "10M" in cigar_str + assert "2D" in cigar_str + + def test_mismatches_ignored(self, sam, tmp_path): + lines = [ + "Seq\tPos\tType\tLen\tRef\tQuery", + "refA_100_aligned_0_F_0_500_0\t10\tmis\t1\tA\tG", + ] + ep = self._write_error_profile(tmp_path, lines) + cigars = sam.get_cigars_nanosim(ep) + assert len(cigars) == 0 # mismatches don't produce entries + + def test_multiple_errors_same_read(self, sam, tmp_path): + lines = [ + "Seq\tPos\tType\tLen\tRef\tQuery", + "refA_100_aligned_0_F_0_500_0\t10\tins\t2\tAC\tACGG", + "refA_100_aligned_0_F_0_500_0\t20\tdel\t3\tACG\t---", + ] + ep = self._write_error_profile(tmp_path, lines) + cigars = sam.get_cigars_nanosim(ep) + cigar_str, _ = cigars["refA-100-0"] + assert "I" in cigar_str + assert "D" in cigar_str + + +# --------------------------------------------------------------------------- +# read_tsv_to_dict +# --------------------------------------------------------------------------- +class TestReadTsvToDict: + def test_basic_parsing(self, sam, tmp_path): + f = tmp_path / "data.tsv" + f.write_text("key1\tval_a\tval_b\nkey2\tval_c\tval_d\n") + result = sam.read_tsv_to_dict(str(f)) + assert result["key1"] == ["val_a", "val_b"] + assert result["key2"] == ["val_c", "val_d"] + + def test_underscore_replaced_with_dash(self, sam, tmp_path): + f = tmp_path / "data.tsv" + f.write_text("my_key_name\tval1\tval2\n") + result = sam.read_tsv_to_dict(str(f)) + assert "my-key-name" in result + + def test_empty_file(self, sam, tmp_path): + f = tmp_path / "empty.tsv" + f.write_text("") + result = sam.read_tsv_to_dict(str(f)) + assert result == {} diff --git a/tests/test_wgsim_to_sam.py b/tests/test_wgsim_to_sam.py new file mode 100644 index 0000000..984db57 --- /dev/null +++ b/tests/test_wgsim_to_sam.py @@ -0,0 +1,42 @@ +"""Tests for pipelines/metagenomic/scripts/wgsim_to_sam.py — CIGAR generation.""" + +import sys +import os +import pytest + +sys.path.insert(0, os.path.join(os.path.dirname(__file__), "..", "pipelines", "metagenomic", "scripts")) +from wgsim_to_sam import get_cigar + + +class TestGetCigar: + def test_all_matches(self): + assert get_cigar("ACGT", "ACGT") == "4M" + + def test_mismatch_counted_as_match(self): + # mismatches are represented by M in CIGAR + assert get_cigar("ACGT", "AGGT") == "4M" + + def test_deletion_in_read(self): + assert get_cigar("ACGT", "A-GT") == "1M1D2M" + + def test_insertion_in_read(self): + assert get_cigar("A-GT", "ACGT") == "1M1I2M" + + def test_consecutive_deletions(self): + assert get_cigar("ACGT", "A--T") == "1M2D1M" + + def test_consecutive_insertions(self): + assert get_cigar("A--T", "ACGT") == "1M2I1M" + + def test_single_base(self): + assert get_cigar("A", "A") == "1M" + + def test_complex_alignment(self): + # 2M + 1D + 1I + 2M (zip pairs: AC match, G/- is D, -/A is I, TT match) + assert get_cigar("ACG-TT", "AC-ATT") == "2M1D1I2M" + + def test_all_deletions(self): + assert get_cigar("AAAA", "----") == "4D" + + def test_all_insertions(self): + assert get_cigar("----", "AAAA") == "4I"