Spaces:
Running
Running
| """Reading sequences and count tables (no model, no cmsearch needed).""" | |
| import gzip | |
| import sys | |
| from pathlib import Path | |
| import pytest | |
| sys.path.insert(0, str(Path(__file__).resolve().parent.parent)) | |
| from webapp.inputs import ( # noqa: E402 | |
| MAX_SEQUENCES, | |
| InputError, | |
| counts_from_size_annotations, | |
| decode_upload, | |
| match_count_table, | |
| parse_sequences, | |
| rna_to_dna, | |
| strip_size, | |
| ) | |
| A = "TACGTAGGGGGCAAGCGTTATCCGGATTTACTGGGTGTAAAGGG" | |
| B = "TACGGAGGATCCGAGCGTTATCCGGATTTATTGGGTTTAAAGGG" | |
| C = "TACGTAGGTGGCAAGCGTTGTCCGGATTTACTGGGCGTAAAGCG" | |
| def _seqs(*pairs): | |
| text = "".join(f">{i}\n{s}\n" for i, s in pairs) | |
| return parse_sequences(text) | |
| # ---------------------------------------------------------------- sequences | |
| def test_fasta_ids_are_the_first_word_as_on_the_cli(): | |
| seqs = parse_sequences(">ASV1 some description\nACGT\nACGT\n>ASV2\nTTTT\n") | |
| assert [(s[0], s[1]) for s in seqs] == [("ASV1", "ACGTACGT"), ("ASV2", "TTTT")] | |
| def test_plain_sequences_get_numbered_names(): | |
| seqs = parse_sequences("ACGT\n\nTTTT\n") | |
| assert [s[0] for s in seqs] == ["seq1", "seq3"] | |
| def test_gzip_upload_is_inflated(): | |
| raw = b">a\nACGT\n" | |
| assert decode_upload(gzip.compress(raw)) == raw.decode() | |
| def test_other_compression_is_refused_clearly(): | |
| with pytest.raises(InputError, match="gzip"): | |
| decode_upload(b"PK\x03\x04rest") | |
| def test_fastq_is_refused_with_advice(): | |
| with pytest.raises(InputError, match="FASTQ"): | |
| parse_sequences("@read1\nACGT\n+\nIIII\n") | |
| def test_header_without_sequence_is_reported(): | |
| with pytest.raises(InputError, match="no sequence"): | |
| parse_sequences(">a\nACGT\n>b\n") | |
| def test_sequence_limit(): | |
| text = "".join(f">s{i}\nA\n" for i in range(MAX_SEQUENCES + 1)) | |
| with pytest.raises(InputError, match="20,000"): | |
| parse_sequences(text) | |
| def test_size_annotations_become_a_one_sample_table(): | |
| seqs = _seqs(("ASV1;size=120", A), ("ASV2;size=7;", B)) | |
| t = counts_from_size_annotations(seqs) | |
| assert t.samples == ["size"] and t.counts == {0: [120.0], 1: [7.0]} | |
| assert strip_size("ASV2;size=7;") == "ASV2" | |
| def test_size_annotations_need_every_header(): | |
| assert counts_from_size_annotations(_seqs(("ASV1;size=3", A), ("ASV2", B))) is None | |
| # ---------------------------------------------------------------- count tables | |
| SEQS = [("ASV1", A), ("ASV2", B), ("ASV3", C)] | |
| def test_features_as_rows_tsv_with_biom_comment(): | |
| text = "# Constructed from biom file\n#OTU ID\tS1\tS2\nASV1\t10\t0\nASV2\t5\t5\nASV3\t0\t1\n" | |
| t = match_count_table(text, _seqs(*SEQS)) | |
| assert t.orientation == "features-as-rows" and t.keyed_by == "id" | |
| assert t.samples == ["S1", "S2"] | |
| assert t.counts == {0: [10, 0], 1: [5, 5], 2: [0, 1]} | |
| assert t.warnings == [] | |
| def test_dada2_seqtab_samples_as_rows_keyed_by_sequence_csv(): | |
| # write.csv(seqtab): empty corner cell, sequences as column names. | |
| text = f'"","{A}","{B}","{C}"\n"S1",3,4,0\n"S2",1,0,9\n' | |
| t = match_count_table(text, _seqs(*SEQS), "seqtab.csv") | |
| assert t.orientation == "features-as-columns" and t.keyed_by == "sequence" | |
| assert t.samples == ["S1", "S2"] | |
| assert t.counts[0] == [3, 1] and t.counts[2] == [0, 9] | |
| def test_r_write_table_header_one_short(): | |
| text = "S1\tS2\nASV1\t1\t2\nASV2\t3\t4\nASV3\t0\t0\n" | |
| t = match_count_table(text, _seqs(*SEQS)) | |
| assert t.samples == ["S1", "S2"] and t.counts[1] == [3, 4] | |
| def test_size_suffix_in_fasta_still_matches_plain_ids(): | |
| t = match_count_table("id\tS1\nASV1\t4\n", _seqs(("ASV1;size=4", A))) | |
| assert t.counts == {0: [4.0]} | |
| def test_taxonomy_column_is_dropped_with_a_note(): | |
| text = "id\tS1\ttaxonomy\nASV1\t4\tk__Bacteria\nASV2\t1\tk__Bacteria\n" | |
| t = match_count_table(text, _seqs(*SEQS)) | |
| assert t.samples == ["S1"] | |
| assert any("taxonomy" in w for w in t.warnings) | |
| def test_no_matching_ids_is_a_clear_error(): | |
| with pytest.raises(InputError, match="None of the feature IDs"): | |
| match_count_table("id\tS1\nOTU_9\t4\n", _seqs(*SEQS)) | |
| def test_partial_match_warns_with_share_of_counts(): | |
| text = "id\tS1\nASV1\t30\nOTU_X\t10\n" | |
| t = match_count_table(text, _seqs(*SEQS)) | |
| assert t.counts == {0: [30.0]} | |
| joined = " ".join(t.warnings) | |
| assert "25.0% of counts" in joined | |
| assert "2 submitted sequence(s) have no counts" in joined | |
| def test_non_numeric_cell_names_the_row(): | |
| with pytest.raises(InputError, match="row 2"): | |
| match_count_table("id\tS1\nASV1\t3\nASV2\tlots\n", _seqs(*SEQS)) | |
| def test_negative_counts_refused(): | |
| with pytest.raises(InputError, match="negative"): | |
| match_count_table("id\tS1\nASV1\t-3\n", _seqs(*SEQS)) | |
| def test_duplicate_sample_names_refused(): | |
| with pytest.raises(InputError, match="not unique"): | |
| match_count_table("id\tS1\tS1\nASV1\t3\t4\n", _seqs(*SEQS)) | |
| def test_example_files_match_each_other(): | |
| ex = Path(__file__).resolve().parent / "static" / "example" | |
| seqs = parse_sequences((ex / "example_v4.fasta").read_text()) | |
| t = match_count_table((ex / "example_v4_counts.tsv").read_text(), seqs, "example_v4_counts.tsv") | |
| assert len(seqs) == len(t.counts) == 92 | |
| assert len(t.samples) == 6 and t.warnings == [] | |
| # ---------------------------------------------------------------- review fixes | |
| def _fasta(n): | |
| return "".join(f">s{i}\n{A}\n" for i in range(n)).encode() | |
| def test_concatenated_gzip_members_are_all_read(): | |
| text = _fasta(300) | |
| half = text.index(b">s150") | |
| data = gzip.compress(text[:half]) + gzip.compress(text[half:]) | |
| assert len(parse_sequences(decode_upload(data))) == 300 | |
| def test_bgzip_style_blocks_split_mid_record_are_all_read(): | |
| # bgzip writes many members with boundaries anywhere, including inside a header. | |
| text = _fasta(501) | |
| data = b"".join(gzip.compress(text[i:i + 997]) for i in range(0, len(text), 997)) | |
| assert len(parse_sequences(decode_upload(data))) == 501 | |
| def test_truncated_gzip_is_refused(): | |
| data = gzip.compress(_fasta(50)) | |
| with pytest.raises(InputError, match="truncated"): | |
| decode_upload(data[:-12]) | |
| def test_gzip_size_cap_counts_all_members(monkeypatch): | |
| import webapp.inputs as inputs | |
| monkeypatch.setattr(inputs, "MAX_TEXT_BYTES", 1000) | |
| data = gzip.compress(b"A" * 600) + gzip.compress(b"A" * 600) | |
| with pytest.raises(InputError, match="once decompressed"): | |
| decode_upload(data) | |
| def test_rna_is_read_as_dna(): | |
| seqs = parse_sequences(">r\nUACGUAGGGG\n>d\nTACGTAGGGG\n>l\nuacg\n") | |
| out, n = rna_to_dna(seqs) | |
| assert n == 2 | |
| assert [s[1] for s in out] == ["TACGTAGGGG", "TACGTAGGGG", "tacg"] | |
| def test_non_finite_counts_are_refused_with_the_column(value): | |
| with pytest.raises(InputError, match="'S1'.*not a count"): | |
| match_count_table(f"id\tS1\nASV1\t{value}\nASV2\t3\n", _seqs(("ASV1", A), ("ASV2", B))) | |
| def test_non_finite_counts_refused_when_features_are_columns(): | |
| with pytest.raises(InputError, match="not a count"): | |
| match_count_table(f"sample\tASV1\tASV2\nS1\tinf\t3\n", _seqs(("ASV1", A), ("ASV2", B))) | |
| def test_duplicate_feature_ids_are_summed_with_a_warning(): | |
| t = match_count_table("id\tS1\nA1\t1\nA1\t5\n", _seqs(("A1", A))) | |
| assert t.counts == {0: [6.0]} | |
| assert any("more than once in the count table" in w and "'A1'" in w for w in t.warnings) | |
| # ---------------------------------------------------------------- alphabet | |
| def test_iupac_codes_gaps_and_case_are_accepted(): | |
| seqs = parse_sequences(">a\nACGTURYSWKMBDHVN\n>b\nacgt-..-nnry\n>c gapped\nAC GT\tAC\n") | |
| assert [s[0] for s in seqs] == ["a", "b", "c"] | |
| def test_fasta_comment_lines_are_skipped_as_the_package_does(): | |
| seqs = parse_sequences(">a\nACGT\n;old-style comment\nACGT\n") | |
| assert [(s[0], s[1]) for s in seqs] == [("a", "ACGTACGT")] | |
| def test_protein_is_refused_naming_the_record_and_character(): | |
| with pytest.raises(InputError) as e: | |
| parse_sequences(">ok\nACGT\n>p1\nMKVLAAGIV\n>p2\nMSTNPK\n") | |
| msg = str(e.value) | |
| assert "2 sequence(s)" in msg and "'p1'" in msg and "'p2'" in msg | |
| assert "first is 'p1': 'L' at position 4" in msg | |
| assert "Nothing was classified" in msg | |
| def test_a_bad_character_late_in_a_record_is_reported_with_its_position(): | |
| with pytest.raises(InputError, match=r"'b': '\*' at position 5"): | |
| parse_sequences(f">a\n{A}\n>b\nACGT*\n") | |
| def test_plain_text_that_is_not_sequence_is_refused(): | |
| with pytest.raises(InputError, match=r"first is 'seq1': 'e' at position 2"): | |
| parse_sequences("hello world\nthis is not a sequence file\n") | |
| def test_alphabet_error_does_not_guess_what_the_file_holds(): | |
| """The same words for a text file and a protein FASTA: the message states | |
| what is accepted and does not call the input protein.""" | |
| for text in ("hello world\nthis is not a sequence file\n", ">p1\nMKVLAAGIV\n"): | |
| with pytest.raises(InputError) as e: | |
| parse_sequences(text) | |
| msg = str(e.value) | |
| assert "protein" not in msg.lower(), msg | |
| assert "nucleotide sequences (DNA or RNA) only" in msg | |
| def test_mostly_valid_file_is_refused_whole_and_lists_five_records(): | |
| good = "".join(f">g{i}\n{A}\n" for i in range(100)) | |
| bad = "".join(f">bad{i}\nACGTX\n" for i in range(8)) | |
| with pytest.raises(InputError) as e: | |
| parse_sequences(good + bad) | |
| msg = str(e.value) | |
| assert msg.startswith("8 sequence(s) contain characters") | |
| assert all(f"'bad{i}'" in msg for i in range(5)) and "'bad5'" not in msg | |
| assert "and 3 more" in msg | |
| def test_alphabet_error_is_checked_after_the_fastq_and_empty_checks(): | |
| with pytest.raises(InputError, match="FASTQ"): | |
| parse_sequences("@read1\nACGT\n+\nIIII\n") | |
| with pytest.raises(InputError, match="no sequence"): | |
| parse_sequences(">a\nMKV\n>b\n") | |