import numbers
import dataclasses
from ete4 import Tree
from .genome import ExtantGenome, AncestralGenome, Genome, EvolutionaryConceptError
from .iham import IHAM
import abc
[docs]
class AbstractGene(metaclass=abc.ABCMeta):
"""
AbstractGene is an abstract class representing extant or ancestral genes. An AbstractGene is defined by an unique
id, the genome it belongs to and its parent AbstractGene.
Attributes:
| parent (:obj:`pyham.abstractgene.HOG`): Direct parent HOG this gene is descendant to. For and only for top level HOG, parent
attribute is set to None.
| genome (:obj:`pyham.genome.Genome`): Related Genome object.
| arose_by_duplication: pyham.abstractgene.DuplicationNode if this AbstractGene arose by a duplication from its parent otherwise False.
"""
def __init__(self, arose_by_duplication=False, **kwargs):
"""
Attributes:
| arose_by_duplication (:obj:`Bool`): True if this AbstractGene arose by a duplication from its parent.
| ** kwargs: dictionary of attribute and value required to create the: obj:`pyham.abstractgene.AbstractGene`.
"""
self.parent = None
self.genome = None
self._properties = {}
self.arose_by_duplication = arose_by_duplication
[docs]
@abc.abstractmethod
def set_genome(self, genome):
"""Set the genome attribute using the given :obj:`Genome`."""
pass
[docs]
@abc.abstractmethod
def is_singleton(self):
""" Return true if singletons otherwise false"""
pass
[docs]
def search_ancestor_hog_in_ancestral_genome(self, ancestral_genome):
"""
Get the :obj:`HOG` related to the query :obj:`AncestralGenome` from which this :obj:`AbstractGene` is
descendant otherwise return None. In addition, it returns a boolean whether a duplication occurs in between the two
:obj:`AbstractGene`.
Args:
ancestral_genome (:obj:`pyham.genome.AncestralGenome`): Ancestral genome of interest.
Returns:
return :obj:`pyham.abstractgene.HOG` or None and :obj:`bool`
"""
found = None
current_hog = self
paralog = current_hog.arose_by_duplication
# Until we reach the top level HOG or the ancestor HOG is founded
while current_hog.parent is not None and found is None:
current_hog = current_hog.parent
if current_hog.genome == ancestral_genome:
return current_hog, paralog
if current_hog.arose_by_duplication != False:
paralog = True
return found, paralog
[docs]
def get_top_level_hog(self):
"""
Get the related top level :obj:`pyham.abstractgene.HOG`. If :obj:`pyham.abstractgene.HOG` return self, the :obj:`pyham.abstractgene.HOG` is already top level. If :obj:`pyham.abstractgene.Gene` return self, the :obj:`pyham.abstractgene.Gene` is a singleton.
Returns:
return :obj:`pyham.abstractgene.HOG`.
"""
current_hog = self
while current_hog.parent is not None:
current_hog = current_hog.parent
return current_hog
[docs]
def get_at_level(self, genome):
"""
Get for the query HOG, the ancestor or children AbstractGene in the genome of interest
Args:
genome (:obj:`pyham.genome.Genome`): Genome of interest.
Returns:
list of :obj:`pyham.abstractgene.AbstractGene`.
Raises:
TypeError: if genome is not a Genome.
KeyError: If no HOG matches the query genome or it return itself.
"""
if not isinstance(genome, Genome):
raise TypeError("expect subclass obj of '{}', got {}"
.format(Genome.__name__,
type(genome).__name__))
tl = self.get_top_level_hog()
def pre(node, elem):
if node.genome == genome:
elem.append(node)
return elem
def leaf(node, leaf, elem):
if leaf.genome == genome:
elem.append(leaf)
return elem
rlist = tl.visit([], function_extant_gene=leaf, function_prefix=pre)
if not rlist:
raise KeyError("Level {} not found within this HOG".format(genome))
for e in rlist:
if e == self:
raise KeyError("get level cannot return itself".format())
return rlist
def add_property(self, name, value):
self._properties[name] = value
def __getitem__(self, item):
return self._properties[item]
[docs]
class HOG(AbstractGene):
"""
HOG class represents ancestral genes. An AbstractGene is defined the genome it belongs to and its child/parent
AbstractGene. Only top level HOG doesn't have parent attribute (set to None).
Attributes:
| hog_id (:obj:`str`): hog id. Defaults is None.
| children (:obj:`list` of :obj:`pyham.abstractgene.AbstractGene`): A list of direct descendants AbstractGene.
| hogvis (:obj:`pyham.IHAM`): :obj:`pyham.IHAM` object of this HOG.
| duplications (:obj:`list` of :obj:`pyham.abstractgene.DuplicationNode`): list of all duplication node child of this HOG.
"""
def __init__(self, id=None, **kwargs):
super(HOG, self).__init__(**kwargs)
self.hog_id = id
if id is None:
self.hog_id = kwargs.get('og', id) # If we have the gene named in og tag, use this.
self.og = kwargs.get('og')
self.taxon_id = int(kwargs['taxonId']) if 'taxonId' in kwargs else None
self.children = []
self.hogvis = None
self.duplications = []
[docs]
def add_child(self, child_to_add):
"""
Add in the children attribute a given :obj:`pyham.abstractgene.AbstractGene` and update parent attribute in child AbstractGene.
Attributes:
child_to_add (:obj:`pyham.abstractgene.AbstractGene`): add_child to add.
"""
if not isinstance(child_to_add, AbstractGene):
raise TypeError("expect subclass obj of '{}', got {}"
.format(AbstractGene.__name__,
type(child_to_add).__name__))
if self == child_to_add:
raise EvolutionaryConceptError('Cannot be a child of itself')
self.children.append(child_to_add)
child_to_add.parent = self
[docs]
def remove_child(self, child_to_remove):
"""
Remove from the children attribute a given :obj:`AbstractGene` and update parent attribute in child
AbstractGene.
Attributes:
child_to_add (:obj:`AbstractGene`): add_child to add.
"""
if not isinstance(child_to_remove, AbstractGene):
raise TypeError("expect subclass obj of '{}', got {}"
.format(AbstractGene.__name__,
type(child_to_remove).__name__))
if child_to_remove in self.children:
self.children.remove(child_to_remove)
child_to_remove.parent = None
else:
raise ValueError("element not found in the hog children")
[docs]
def score(self, score_id, value=None):
"""
Get or set a score attribute for this HOG.
If no value is passed, the stored value is returned or a KeyError is raised if this score has not been stored.
Attributes:
score_id (:obj:`ExtantGenome`): score to be set or returned.
value (optional): A numeric value for the score of this HOG.
score_id (:obj:`ExtantGenome`): score to be set or returned.
Raises:
KeyError: if score_id is unknown and accessed.
"""
if value is None:
try:
return self.scores[score_id]
except (AttributeError, KeyError):
raise KeyError("Score '{}' not set".format(score_id))
if not isinstance(value, numbers.Number):
raise ValueError("value must be numeric")
scores = getattr(self, 'scores', {})
scores[score_id] = value
self.scores = scores
[docs]
def set_genome(self, genome):
"""
This method set :obj:`pyham.genome.AncestralGenome` given as parameter as genome attribute.
Attributes:
| genome (:obj:`pyham.genome.AncestralGenome`): Ancestral genome to be set as genome.
"""
# Check that genome is an :obj:`AncestralGenome`
if not isinstance(genome, AncestralGenome):
raise TypeError("expect subclass obj of '{}', got {}"
.format(AncestralGenome.__name__,
type(genome).__name__))
if self.genome is not None and genome != self.genome:
raise EvolutionaryConceptError("HOG can only be mapped to one ancestral genome")
self.genome = genome
[docs]
def visit(self, elem, function_extant_gene=None, function_postfix=None, function_prefix=None):
"""
Recursive function to traverse the HOG nested children hierarchy.
An arbitrary object "elem" is return at each visit() and is pass/returned from each recursive state function
calls.
There is tree recursive state functions:
- function_extant_gene: Function called at each ExtantGene child. Its take as argument the child ExtantGene
and the "elem" and return processed "elem".
- function_postfix: Function called after each HOG child visit(). Its take as argument the child HOG and
the "elem" and return processed "elem".
- function_prefix: First function called that take as argument the current HOG and the "elem" and return
processed "elem".
Attributes:
| elem (arbitrary object): object that is pass through visit recursive and recursive state functions calls.
| function_extant_gene (callback function): Callback function for ExtantGenes (leaves).
| function_postfix (callback function): Callback function for after child visit current HOG.
| function_prefix (callback function): Callback function for current HOG.
Returns:
return arbitrary
"""
if function_prefix is not None:
elem = function_prefix(self, elem)
for child in self.children:
if isinstance(child, Gene):
if function_extant_gene is not None:
elem = function_extant_gene(self, child, elem)
else:
child.visit(elem, function_extant_gene=function_extant_gene, function_postfix=function_postfix,
function_prefix=function_prefix)
if function_postfix is not None:
elem = function_postfix(self, child, elem)
return elem
# TODO: UT
[docs]
def get_all_descendant_genes(self):
"""
Get all :obj:`Gene` present in this HOG.
Returns:
list of :obj:`pyham.abstractgene.Gene`
"""
def append_leaf(current, child, list):
list.append(child)
return list
return self.visit([], function_extant_gene=append_leaf)
[docs]
def get_all_descendant_genes_clustered_by_species(self):
"""
Get all :obj:`pyham.abstractgene.Gene` present in this HOG clustered by species :obj:`pyham.genomeExtantGenome`.
Returns:
Dictionary of :obj:`pyham.genome.ExtantGenome` map their list of :obj:`pyham.abstractgene.Gene`.
"""
def append_add_lead_with_genome(current, child, dict):
dict.setdefault(child.genome, []).append(child)
return dict
return self.visit({}, function_extant_gene=append_add_lead_with_genome)
[docs]
def get_all_descendant_hogs(self): # TODO: self is also returned with it
"""
Get all :obj:`pyham.abstractgene.HOG` present in this HOG hierarchy.
Returns:
Dictionary of :obj:`pyham.genomeExtantGenome` map their list of :obj:`pyham.abstractgene.Gene`.
"""
def append_child(current, list):
list.append(current)
return list
return self.visit([], function_prefix=append_child)
[docs]
def get_all_descendant_hog_levels(self): # TODO: self is also returned with it
"""
Get all :obj:`pyham.genome.Genome` present in this HOG hierarchy.
Returns:
List of :obj:`pyham.genome.Genome`.
"""
def append_current_genome(current,list):
list.append(current.genome)
return list
return self.visit([], function_prefix=append_current_genome)
[docs]
def get_number_losses(self):
"""
Compute the number of implied gene losses within this HOG's hierarchy, following the
Dollo parsimony principle: once a gene lineage is lost it cannot be regained, so a
whole clade lacking any representative under this HOG is explained by a single loss
event placed on the branch leading to that clade, rather than one loss per missing
descendant leaf.
For every ancestral level in this HOG's subtree (starting at self), each immediate
child taxon in the species tree that has no matching child :obj:`HOG`/:obj:`Gene` here
contributes one loss; that clade is not explored further since it is already fully
accounted for by this single event.
Returns:
:obj:`int`: number of implied Dollo loss events in this HOG's subtree.
"""
present_children_taxa = {child.genome.taxon for child in self.children}
losses = sum(1 for taxon in self.genome.taxon.children if taxon not in present_children_taxa)
for child in self.children:
if isinstance(child, HOG):
losses += child.get_number_losses()
return losses
[docs]
def get_hog_vis(self, newick_str):
""" Lazy getter of the :obj:`pyham.abstractgene.HOG` IHAM.
Args:
newick_str (:obj:`str`): newick species tree used by the create_iHam.
Returns:
:obj:`pyham.IHAM` of this HOG.
"""
if self.hogvis is None:
self.hogvis = IHAM(newick_str, self)
return self.hogvis
# TODO: end UT
[docs]
def is_singleton(self):
return False
[docs]
def to_ete_tree(self):
"""
Convert this :obj:`HOG` and its nested sub-HOGs into an :obj:`ete4.Tree`.
Each nested :obj:`HOG` becomes an internal node named after its hog_id
(also stored as the "hog_id" prop), and each terminal :obj:`Gene` becomes
a leaf named after its external gene id (geneId, protId or transcriptId,
in that order of preference, falling back to its internal unique_id),
with all available external ids attached as leaf props.
Returns:
:obj:`ete4.Tree`: tree representation of this HOG hierarchy, rooted at self.
"""
def leaf_name(gene):
xref = gene.get_dict_xref()
for key in ("geneId", "protId", "transcriptId", "id"):
if key in xref:
return str(xref[key])
return str(gene.unique_id)
def build(hog):
node = Tree({"name": hog.hog_id})
node.add_prop("hog_id", hog.hog_id)
node.add_prop("taxon", hog.genome.name)
node.add_prop("duplicate", bool(hog.arose_by_duplication))
if hog.arose_by_duplication:
# identifies which sibling nodes share the same duplication event
# (hog.arose_by_duplication is a shared DuplicationNode instance);
# a plain object id is enough since this only needs to be stable
# within one export, not meaningful across separate runs.
node.add_prop("duplication_id", id(hog.arose_by_duplication))
for child in hog.children:
if isinstance(child, Gene):
leaf = node.add_child(name=leaf_name(child))
for key, value in child.get_dict_xref().items():
leaf.add_prop(key, value)
leaf.add_prop("taxon", child.genome.name)
else:
node.add_child(build(child))
return node
return build(self)
[docs]
def find_by_id(self, hog_id: str, taxid: int = None):
"""
Recursive function to find a HOG by its id and optionally its taxid.
Args:
| hog_id (:obj:`str`): HOG id to be found.
| taxid (optional, :obj:`int`): taxid of the genome of the HOG to be found.
Returns:
:obj:`pyham.abstractgene.HOG` if found otherwise None.
"""
cur_hog_without_taxid = self.hog_id.split('_')[0]
if cur_hog_without_taxid == hog_id and (taxid is None or self.taxon_id == taxid or self.hog_id == f"{hog_id}_{taxid}"):
return self
elif not hog_id.startswith(cur_hog_without_taxid):
return None
for child in getattr(self, 'children', []):
if isinstance(child, HOG):
found = child.find_by_id(hog_id, taxid)
if found is not None:
return found
return None
def __repr__(self):
ids = []
if self.hog_id is not None:
ids.append("id={}".format(self.hog_id))
if self.og is not None and self.og != self.hog_id:
ids.append("og={}".format(self.og))
if len(ids) == 0:
n = self.parent
while n.hog_id is None and n.parent is not None:
n = n.parent
ids.append("id={}-{}".format(n.hog_id, self.__hash__()))
try:
ids.append("level={}".format(self.genome.name))
except AttributeError:
pass
return "<{}({})>".format(self.__class__.__name__, ";".join(ids))
[docs]
class Gene(AbstractGene):
"""Gene class for extant genes (inherit from AbstractGene).
Attributes:
| unique_id internal unique id
| gene_id id used to mapped to external ids.
| prot_id id used to mapped to external ids.
| transcript_id id used to mapped to external ids.
"""
def __init__(self, id, geneId=None, protId=None, transcriptId=None, **kwargs):
super(Gene, self).__init__(**kwargs)
self.unique_id = id
self.gene_id = geneId
self.prot_id = protId
self.transcript_id = transcriptId
[docs]
def set_LOFT(self, loft_id):
"""Set the LOFT identifier for the gene.
In case the ID has already been set, a ValueError is raised. It is assumed
that a gene can only have one LOFT identifier."""
if hasattr(self, 'hog_id'):
raise ValueError("Gene {} has already a LOFT identifier assigned ({}).".format(self, self.hog_id))
self.hog_id = loft_id
[docs]
def set_genome(self, genome):
"""
This method set :obj:`pyham.genome.ExtantGenome` given as parameter as genome attribute.
Attributes:
genome (:obj:`pyham.genome.ExtantGenome`): Extant genome to be set as genome.
Raises:
TypeError: If genome is not :obj:`pyham.genome.ExtantGenome`.
EvolutionaryConceptError: If we try to add more than one genome.
"""
# Check that genome is an :obj:`ExtantGenome`
if not isinstance(genome, ExtantGenome):
raise TypeError("expect subclass obj of '{}', got {}"
.format(ExtantGenome.__name__,
type(genome).__name__))
if self.genome is not None and genome != self.genome:
raise EvolutionaryConceptError("Gene can only belong to one genome")
self.genome = genome
[docs]
def get_dict_xref(self): # <-- TODO: UT
""" Get the dictionary of cross references.
Returns:
Dictionary of external id tag map with their value.
"""
xref = {}
if self.unique_id is not None:
xref["id"] = self.unique_id
if self.gene_id is not None:
xref["geneId"] = self.gene_id
if self.prot_id is not None:
xref["protId"] = self.prot_id
if self.transcript_id is not None:
xref["transcriptId"] = self.transcript_id
return xref
[docs]
def is_singleton(self):
if self.parent is None:
return True
else:
return False
def __repr__(self):
return "{}({})".format(self.__class__.__name__, self.unique_id)
def short_level_name(name, max_len=60):
"""Abbreviate a level name synthesized from its children's names (``A/B/C/...``) for messages."""
if name is None or len(name) <= max_len or "/" not in name:
return name
parts = name.split("/")
return f"{parts[0]}/.../{parts[-1]} ({len(parts)} species)"
class SameLevelHOGError(Exception):
pass
@dataclasses.dataclass
class TaxonomicConflict:
"""One HOG whose claimed taxonomic level conflicts with the given species tree.
Attributes:
| family_id (:obj:`str`): top-level HOG id this conflict was found in.
| hog_id (:obj:`str`): the HOG whose TaxRange/taxid claim is in question.
| claim_property (:obj:`str`): "TaxRange" or "taxid".
| claim_value (:obj:`str`): the raw property value that was claimed.
| resolved_level (:obj:`str`): name of the taxon node actually used. Equal to
claim_value at present (this only fires once the claim itself resolved to a
node), kept as a separate field for forward-compatibility rather than as a
redundant duplicate.
| kind (:obj:`str`): "same_rank" (a child collapses onto the same level as
its own parent), "disjoint" (a child sits in an unrelated branch), or
"inverted" (the claimed level is itself a descendant of the child).
| offending_members (:obj:`list`): list of (label, species, detail) tuples
describing each child that doesn't fit under the claimed level.
| sibling_members (:obj:`list`): list of (label, species) tuples describing
this HOG's other children, for context on what else is present at this level.
"""
family_id: str
hog_id: str
claim_property: str
claim_value: str
resolved_level: str
kind: str
offending_members: list
sibling_members: list
class TaxonomicConflictError(Exception):
"""Raised when one or more HOGs' claimed taxonomic level conflicts with the
given species tree. Carries the full list of :obj:`TaxonomicConflict` records
found (possibly just one, in fail_fast mode) in the `conflicts` attribute."""
def __init__(self, message: str, conflicts: list):
super().__init__(message)
self.conflicts = conflicts
class DuplicationNode(object):
"""
This object link together all abstract genes that emerges from the same duplication event. Its composed of a set of genes
all descendant from this duplication event and their direct related parent hog. In addition it contained the parent of the MRCA of all
children genes (otherwise MRCA return genome father duplication).
Attributes:
| parent (:obj:`pyham.abstractgene.HOG`): Direct parent HOG this DuplicationNode.
| children: List of all duplicated genes
| MRCA: Parent of the MRCA of all children gene that contained the parent HOG.
"""
def __init__(self, parser, id=None):
self.parser = parser
self.MRCA = None
self.children = []
self.parent = None
self.id = id
def set_parent(self, hog):
"""
Args:
| hog (:obj:`pyham.abstractegene.HOG`): parent :obj:`pyham.abstractegene.HOG` of this duplication.
"""
if not isinstance(hog, HOG):
raise TypeError("expect subclass obj of '{}', got {}"
.format(HOG.__name__,
type(hog).__name__))
hog.duplications.append(self)
self.parent = hog
def set_MRCA(self):
"""
Compute the MRCA of all genes genomes.
It is always defined as the level of the genome which is the parent
of the branch on which the duplication event occured. E.g. if the
duplication event on the branch between Mammalia and Euchontoglires,
the MRCA is set to Mammalia.
"""
children_taxa = set([child.genome.taxon for child in self.children])
if len(children_taxa) < 2:
mrca = children_taxa.pop().up
else:
mrca = self.parser.taxonomy.get_mrca_taxnode(*children_taxa)
if mrca.up is None:
copies = ", ".join(
f"{c.hog_id if isinstance(c, HOG) else (c.prot_id or c.gene_id or c.unique_id)} "
f"at \"{short_level_name(c.genome.name)}\"" for c in self.children)
raise EvolutionaryConceptError(
f"The copies of {'duplication ' + str(self.id) if self.id is not None else 'a duplication'} "
f"({copies}) only meet at the root of the species tree (\"{short_level_name(mrca.name)}\"), "
f"so the duplication would predate the tree. The species tree probably lacks levels the "
f"OrthoXML was built with, e.g. because it is pruned.")
mrca = mrca.up
self.MRCA = self.parser.taxonomy.get_genome_from_taxnode(mrca)
def add_child(self, child):
"""
Add child to list of children
"""
if not isinstance(child, AbstractGene):
raise TypeError("expect subclass obj of '{}', got {}"
.format(AbstractGene.__name__,
type(child).__name__))
self.children.append(child)
child.arose_by_duplication = self
def remove_child(self, child_to_remove):
"""
remove child to list of children
"""
if not isinstance(child_to_remove, AbstractGene):
raise TypeError("expect subclass obj of '{}', got {}"
.format(AbstractGene.__name__,
type(child_to_remove).__name__))
if child_to_remove in self.children:
self.children.remove(child_to_remove)
child_to_remove.arose_by_duplication = False
else:
raise ValueError("element not found in the duplication children")