Skip to content

New source: SILVA taxonomy #348

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
wants to merge 18 commits into
base: main
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from 7 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 2 additions & 0 deletions src/pyobo/sources/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -58,6 +58,7 @@
from .selventa import SCHEMGetter, SCOMPGetter, SDISGetter, SFAMGetter
from .sgd import SGDGetter
from .signor import SignorGetter
from .silva import SILVAGetter
from .slm import SLMGetter
from .umls import UMLSGetter, UMLSSTyGetter
from .unimod import UnimodGetter
Expand Down Expand Up @@ -139,6 +140,7 @@
"SDISGetter",
"SFAMGetter",
"SGDGetter",
"SILVAGetter",
"SLMGetter",
"SignorGetter",
"UMLSGetter",
Expand Down
173 changes: 173 additions & 0 deletions src/pyobo/sources/silva.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,173 @@
"""Convert SILVA small subunit (ssu) taxonomy to OBO format."""

import logging
from collections.abc import Iterable

import pandas as pd
from tqdm.auto import tqdm

from pyobo.struct import Obo, Reference, Term, TypeDef, default_reference
from pyobo.struct.typedef import has_taxonomy_rank
from pyobo.utils.path import ensure_path

__all__ = [
"SILVAGetter",
]

PREFIX = "silva.taxon"

#: A mapping from SILVA rank names to TAXRANK references
SILVA_RANK_TO_TAXRANK = {
"domain": Reference(prefix="TAXRANK", identifier="0000037", name="domain"),
"major_clade": Reference(prefix="TAXRANK", identifier="0001004", name="major_clade"),
"superkingdom": Reference(prefix="TAXRANK", identifier="0000022", name="superkingdom"),
"kingdom": Reference(prefix="TAXRANK", identifier="0000017", name="kingdom"),
"subkingdom": Reference(prefix="TAXRANK", identifier="0000029", name="subkingdom"),
"superphylum": Reference(prefix="TAXRANK", identifier="0000027", name="superphylum"),
"phylum": Reference(prefix="TAXRANK", identifier="0000001", name="phylum"),
"subphylum": Reference(prefix="TAXRANK", identifier="0000008", name="subphylum"),
"infraphylum": Reference(prefix="TAXRANK", identifier="0000040", name="infraphylum"),
"superclass": Reference(prefix="TAXRANK", identifier="0000015", name="superclass"),
"class": Reference(prefix="TAXRANK", identifier="0000002", name="class"),
"subclass": Reference(prefix="TAXRANK", identifier="0000007", name="subclass"),
"infraclass": Reference(prefix="TAXRANK", identifier="0000019", name="infraclass"),
"superorder": Reference(prefix="TAXRANK", identifier="0000020", name="superorder"),
"order": Reference(prefix="TAXRANK", identifier="0000003", name="order"),
"suborder": Reference(prefix="TAXRANK", identifier="0000014", name="suborder"),
"superfamily": Reference(prefix="TAXRANK", identifier="0000018", name="superfamily"),
"family": Reference(prefix="TAXRANK", identifier="0000004", name="family"),
"subfamily": Reference(prefix="TAXRANK", identifier="0000024", name="subfamily"),
"genus": Reference(prefix="TAXRANK", identifier="0000005", name="genus"),
}

#: URLs for the SILVA files.
SILVA_TAXONOMY_URL = "https://www.arb-silva.de/fileadmin/silva_databases/current/Exports/taxonomy/tax_slv_ssu_138.2.txt.gz"
SILVA_TAXMAP_URL = "https://www.arb-silva.de/fileadmin/silva_databases/current/Exports/taxonomy/taxmap_slv_ssu_ref_nr_138.2.txt.gz"

logger = logging.getLogger(__name__)
logger.setLevel(logging.WARNING)

TYPEDEF = TypeDef(
reference=default_reference(PREFIX, "fixme", name="fixme"),
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

fixme here

definition="fixme",
is_metadata_tag=False,
)


class SILVAGetter(Obo):
"""An ontology representation of the SILVA taxonomy."""

ontology = bioversions_key = PREFIX
typedefs = [has_taxonomy_rank, TYPEDEF]
idspaces = {
PREFIX: "https://www.arb-silva.de/no_cache/download/archive/current/Exports/taxonomy/",
"ena.embl": "https://www.ebi.ac.uk/ena/browser/view/",
}
root_terms = [
Reference(prefix=PREFIX, identifier="2", name="Archaea"),
Reference(prefix=PREFIX, identifier="3", name="Bacteria"),
Reference(prefix=PREFIX, identifier="4", name="Eukaryota"),
]

def iter_terms(self, force: bool = False) -> Iterable[Term]:
"""Iterate over terms in the SILVA ontology."""
return iter_terms_silva(version=self._version_or_raise, force=force)

Check warning on line 74 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L74

Added line #L74 was not covered by tests


def iter_terms_silva(version: str, force: bool = False) -> Iterable[Term]:
"""Iterate over SILVA terms from the main taxonomy and taxmap files."""
# --- Process the main taxonomy file ---
taxonomy_path = ensure_path(PREFIX, url=SILVA_TAXONOMY_URL, version=version, force=force)
tax_df = pd.read_csv(

Check warning on line 81 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L80-L81

Added lines #L80 - L81 were not covered by tests
taxonomy_path,
sep="\t",
header=None,
names=["taxonomy", "taxon_id", "rank", "ignore", "introduced"],
dtype=str,
)
tax_df.fillna("", inplace=True)

Check warning on line 88 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L88

Added line #L88 was not covered by tests

#: a dictionary that maps the joined taxonomy path (with trailing ";") to taxon_id
tax_path_to_id: dict[str, str] = {}

Check warning on line 91 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L91

Added line #L91 was not covered by tests
#: maps taxon_id to the Term object
terms_by_id = {}

Check warning on line 93 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L93

Added line #L93 was not covered by tests

for idx, row in tqdm(
tax_df.iterrows(),
total=len(tax_df),
desc=f"[{PREFIX}] processing main taxonomy",
unit="row",
):
tax_str = row["taxonomy"].strip()
taxon_id = row["taxon_id"].strip()
rank_raw = row["rank"].strip()
rank = rank_raw.lower()

Check warning on line 104 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L101-L104

Added lines #L101 - L104 were not covered by tests
# Split taxonomy string by ";" and discard empty parts.
parts = [p.strip() for p in tax_str.split(";") if p.strip()]

Check warning on line 106 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L106

Added line #L106 was not covered by tests
if not parts:
logger.warning(f"Row {idx}: empty taxonomy string: {tax_str}")
continue

Check warning on line 109 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L108-L109

Added lines #L108 - L109 were not covered by tests

# The term's name is the last element (e.g. for "Bacteria;Actinomycetota;", name is "Actinomycetota").
name = parts[-1]
term = Term(reference=Reference(prefix=PREFIX, identifier=taxon_id, name=name))

Check warning on line 113 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L112-L113

Added lines #L112 - L113 were not covered by tests
if rank in SILVA_RANK_TO_TAXRANK:
term.annotate_object(has_taxonomy_rank, SILVA_RANK_TO_TAXRANK[rank])

Check warning on line 115 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L115

Added line #L115 was not covered by tests
else:
logger.warning(

Check warning on line 117 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L117

Added line #L117 was not covered by tests
f"Row {idx}: unknown rank '{rank_raw}' for taxonomy: {tax_str} (taxon id: {taxon_id})"
)

# Determine the parent by joining all but the last element.
if len(parts) > 1:
parent_key = ";".join(parts[:-1]) + ";" # e.g. "Bacteria;"
parent_id = tax_path_to_id.get(parent_key)

Check warning on line 124 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L123-L124

Added lines #L123 - L124 were not covered by tests
if parent_id:
term.append_parent(Reference(prefix=PREFIX, identifier=parent_id))
full_key = ";".join(parts) + ";"
tax_path_to_id[full_key] = taxon_id
terms_by_id[taxon_id] = term

Check warning on line 129 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L126-L129

Added lines #L126 - L129 were not covered by tests

# --- Process the taxmap file ---
# This file has a header with columns: primaryAccession, start, stop, path, organism_name, taxid
taxmap_path = ensure_path(PREFIX, url=SILVA_TAXMAP_URL, version=version, force=force)
taxmap_df = pd.read_csv(taxmap_path, sep="\t", dtype=str)
taxmap_df.rename(

Check warning on line 135 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L133-L135

Added lines #L133 - L135 were not covered by tests
columns={
"primaryAccession": "accession",
"organism_name": "organism",
"taxid": "species_taxon_id",
"path": "taxonomy",
},
inplace=True,
)
taxmap_df.fillna("", inplace=True)

Check warning on line 144 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L144

Added line #L144 was not covered by tests

for idx, row in tqdm(
taxmap_df.iterrows(), total=len(taxmap_df), desc=f"[{PREFIX}] processing taxmap", unit="row"
):
accession = row["accession"].strip()
species_taxon_id = row["species_taxon_id"].strip()
organism = row["organism"].strip()

Check warning on line 151 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L149-L151

Added lines #L149 - L151 were not covered by tests
if not accession or not species_taxon_id:
continue

Check warning on line 153 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L153

Added line #L153 was not covered by tests
if species_taxon_id in terms_by_id:
# Create a new term for the ENA accession.
new_term = Term(

Check warning on line 156 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L156

Added line #L156 was not covered by tests
reference=Reference(prefix="ena.embl", identifier=accession, name=organism)
)
# Do NOT annotate the new term with a rank (leave it unranked).
new_term.append_parent(Reference(prefix=PREFIX, identifier=species_taxon_id))
Copy link
Member

@cthoyt cthoyt Feb 13, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Don't ENA terms represent nucleotide sequences derived from experiments? Can they also represent projects?

From what I understand, they aren't actually themselves representing taxa. Therefore this parent/child relationship doesn't make sense.

The hard work of making a PyOBO source is really understanding what is the relationship SILVA means when it mentions its internal taxonomy and ENA sequences. I can't do this hard work for you in detail, but from a high level it seems like the sequence was derived from an individual of the taxonomy.

Then, there's two options:

  1. Find an existing RO relationship that is appropriate for this. Maybe http://purl.obolibrary.org/obo/RO_0001001, even though it's not a perfect ontological fit. Maybe OBI is a better place to look
  2. mint an ad-hoc one yourself within the scope of this file, e.g., like in
    HAS_INTERVENTION = TypeDef(

If you go the second route, make sure that you do a good job describing what the relationship means (in a concise way)

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thank you for your detailed feedback. I completely understand where the hard work lies, and I truly appreciate the guidance you provided. Your suggestions—either reusing an existing RO relationship (like RO_0001001) or minting an ad-hoc one (as in clinicaltrials.py)—are exactly the direction I was hoping for.

I’ll explore those options further. Alternatively, I might start by representing only down to the genus level (as shown in the taxonomy files) until I fully understand the nuances of the lower levels.

Thanks again for steering this work in the right direction!

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think it's good to surface decision like this higher up, and ideally all of pyobo would be biolink compliant. biolink:has_biological_sequence is the right KG relationship to use.

If you are going to use RO you need to use it consistently with how it's intended and not just pick a label that sounds right.

yield new_term

Check warning on line 161 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L160-L161

Added lines #L160 - L161 were not covered by tests
else:
logger.warning(

Check warning on line 163 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L163

Added line #L163 was not covered by tests
f"Row {idx} in taxmap: species_taxon_id {species_taxon_id} not found in main taxonomy"
)

# Yield all terms from the main taxonomy.
for term in terms_by_id.values():
yield term

Check warning on line 169 in src/pyobo/sources/silva.py

View check run for this annotation

Codecov / codecov/patch

src/pyobo/sources/silva.py#L169

Added line #L169 was not covered by tests


if __name__ == "__main__":
SILVAGetter().cli()