From ea72dd5df390caa78ac754de6463520eda4d4ef7 Mon Sep 17 00:00:00 2001 From: Oliver Kohlbacher Date: Tue, 1 Sep 2026 12:35:47 +0200 Subject: [PATCH 1/4] Add ProForma 2.0 peptide exporter (-epepproforma) Exports each peptide path as a ProForma 2.0 string (HUPO-PSI standard) into a TSV of accession, start, end, misscleavages, proforma. Modifications from -fm/-vm are emitted as signed mass-delta tags by default, optionally as named tags via --pep_proforma_mod_names; terminal modifications use ProForma terminal syntax. Placement is edge-based (first residue of the traversed edge's target node), which also covers modifications on residues introduced by VARIANT/MUTAGEN/CONFLICT, whose features carry no reference location. Collapsed-edge Or qualifiers export only modifications shared by all branches, with a warning when they differ. Co-Authored-By: Claude Fable 5 --- README.md | 22 ++++ protgraph/cli.py | 46 +++++++ protgraph/export/exporters.py | 1 + protgraph/export/peptides/pep_proforma.py | 146 ++++++++++++++++++++++ protgraph/protgraph.py | 1 + tests/fixtures/minimal.txt | 16 +++ tests/functional_test.py | 50 ++++++++ 7 files changed, 282 insertions(+) create mode 100644 protgraph/export/peptides/pep_proforma.py create mode 100644 tests/fixtures/minimal.txt diff --git a/README.md b/README.md index 4eb2b8f..01d1e31 100644 --- a/README.md +++ b/README.md @@ -366,6 +366,28 @@ $ cat e_coli_with_selected_features_limited.fasta| grep "^>" | wc -l 648551 ``` +### Exporting peptides as ProForma 2.0 + +Peptides can also be exported in [ProForma 2.0](https://github.com/HUPO-PSI/ProForma) +notation (the HUPO-PSI standard for writing peptidoforms/proteoforms) into a +tab-separated file with the columns accession, start, end, misscleavages and +proforma: + +```shell +$ protgraph -d trypsin --pep_hops 3 -fm "C:57.021464" -vm "M:15.994915" -epepproforma --pep_proforma_out peptides.proforma.tsv examples/e_coli.dat +$ grep -m 2 "\[" peptides.proforma.tsv +P0A884 49 53 1 RC[+57.021464]HLR +P0A884 1 2 0 M[+15.994915]K +``` + +Modifications applied via `-fm`/`-vm` are known to ProtGraph only by their +delta mass, so they are written as mass-delta tags (`C[+57.021464]`), which is +valid ProForma. If you want named tags instead, declare a mapping for a delta +you passed: `--pep_proforma_mod_names "15.994915=UNIMOD:35"` then writes +`M[UNIMOD:35]`. Terminal modifications (`NPEPTERM` etc.) use the ProForma +terminal syntax (`[+42.010565]-PEPTIDE`). + + We can see that the number peptides within this FASTA is significantly reduced (to ~650 000 entries). **NOTE:** for setting an upper length-limit of peptides we need to set `-nm`. In case of not setting this parameter, longer peptides may be exported. Finally, since the FASTA-file still contains for some entries same sequences, we can concatenate these by using a different (and more sophisticated) FASTA-exporter: diff --git a/protgraph/cli.py b/protgraph/cli.py index 63f4d18..7c024cb 100644 --- a/protgraph/cli.py +++ b/protgraph/cli.py @@ -662,6 +662,52 @@ def add_fasta_peptide_export(group): ) +def add_proforma_peptide_export(group): + group.add_argument( + "--export_peptide_proforma", "-epepproforma", default=False, action="store_true", + help="Set this flag to export peptides as ProForma 2.0 strings (PSI standard) " + "into a single tab-separated file (accession, start, end, misscleavages, proforma). " + "Modifications from -fm/-vm are written as mass-delta tags (e.g. 'S[+79.966]')." + ) + group.add_argument( + "--pep_proforma_out", default=os.path.join(os.getcwd(), "peptides.proforma.tsv"), + type=str, + help="Set the output file for the ProForma peptide export. " + "Default: '${pwd}/peptides.proforma.tsv'. NOTE: This will overwrite existing files." + ) + + def _parse_mod_name(input: str): + if "=" not in input: + raise ArgumentTypeError( + "A modification name mapping needs the form '=', " + "e.g. '79.966=UNIMOD:21'. Received: '{}'".format(input) + ) + delta, name = input.split("=", 1) + try: + delta = float(delta) + finite = delta == delta and abs(delta) != float("inf") + except ValueError: + finite = False + if not finite: + raise ArgumentTypeError( + "The DeltaMass of a name mapping is not a finite number. Received: '{}'".format(input) + ) + name = name.strip() + if not name or any(c in name for c in "[]\t\n\r"): + raise ArgumentTypeError( + "The Name of a name mapping is empty or contains brackets/whitespace " + "that would corrupt the output. Received: '{}'".format(input) + ) + return delta, name + group.add_argument( + "--pep_proforma_mod_names", type=_parse_mod_name, action="append", default=None, + help="Optionally map a delta mass used in -fm/-vm to a named ProForma tag, " + "e.g. --pep_proforma_mod_names '79.966=UNIMOD:21' writes 'S[UNIMOD:21]' instead " + "of 'S[+79.966]'. Can be provided multiple times. No mass lookup is performed; " + "the name is taken as declared." + ) + + def add_trie_peptide_export(group): group.add_argument( "--export_peptide_trie", "-epeptrie", default=False, action="store_true", diff --git a/protgraph/export/exporters.py b/protgraph/export/exporters.py index 131410d..da36ed6 100644 --- a/protgraph/export/exporters.py +++ b/protgraph/export/exporters.py @@ -79,6 +79,7 @@ def __available_exporters(self, **kwargs): # Peptide Exporter to local filesystem (no setup required) (kwargs["export_peptide_fasta"], "protgraph.export.peptides.pep_fasta", "PepFasta"), + (kwargs["export_peptide_proforma"], "protgraph.export.peptides.pep_proforma", "PepProForma"), (kwargs["export_peptide_trie"], "protgraph.export.peptides.pep_trie", "PepTrie"), (kwargs["export_peptide_sqlite"], "protgraph.export.peptides.pep_sqlite", "PepSQLite"), ] diff --git a/protgraph/export/peptides/pep_proforma.py b/protgraph/export/peptides/pep_proforma.py new file mode 100644 index 0000000..588d421 --- /dev/null +++ b/protgraph/export/peptides/pep_proforma.py @@ -0,0 +1,146 @@ +from decimal import Decimal +from math import isfinite + +from protgraph.export.peptides.pep_fasta import PepFasta +from protgraph.graph_collapse_edges import Or + + +class PepProForma(PepFasta): + """ + Peptide exporter emitting ProForma 2.0 strings (HUPO-PSI standard notation + for proteoforms and peptidoforms, LeDuc et al., J. Proteome Res. 2022). + + Writes one tab-separated line per exported peptide path: + + accession start end misscleavages proforma + + ProtGraph knows modifications only by their delta mass (-fm/-vm), so + modifications are emitted as mass-delta tags by default ("S[+79.966]"), + which is valid ProForma. A user-declared mapping (--pep_proforma_mod_names + "79.966=UNIMOD:21") replaces a delta by a named tag ("S[UNIMOD:21]"). + Terminal modifications (NPEPTERM etc.) use ProForma terminal syntax + ("[+42.010565]-PEPTIDE"). + + A modification is placed on the residue its edge points at: FIXMOD/VARMOD + features are attached to the in-edges of the (cloned) modified node when + the graph is generated, so the first residue of the traversed edge's + target node is the modified one. This also holds for residues introduced + by VARIANT/MUTAGEN/CONFLICT, whose features carry no usable reference + location, and for nodes merged after annotation. + + NOTE: as with the FASTA exporter, this exports all possible paths, so make + sure the traversal can terminate in forseeable future! + """ + + def start_up(self, **kwargs): + super(PepProForma, self).start_up(**kwargs) + self.output_file = kwargs["pep_proforma_out"] + self.mod_names = dict(kwargs["pep_proforma_mod_names"] or []) + self.warned_or_once = False + + def export_peptides(self, prot_graph, l_path_nodes, l_path_edges, l_peptide, l_miscleavages, queue): + entries = "" + for peptide, nodes, edges, misses in zip(l_peptide, l_path_nodes, l_path_edges, l_miscleavages): + # the first/last node carrying sequence; with terminal modifications + # applied, nodes[1] can be an empty helper node (see annotate_ptms) + inner = [n for n in nodes[1:-1] if prot_graph.vs[n]["aminoacid"]] + acc = self._get_accession_or_isoform(prot_graph.vs[inner[0]]) + start_pos = self._get_position_or_isoform_position(prot_graph.vs[inner[0]]) + end_pos = self._get_position_or_isoform_position(prot_graph.vs[inner[-1]], end=True) + proforma = self._build_proforma(prot_graph, inner, edges, peptide) + entries += "\t".join( + [acc, str(start_pos), str(end_pos), str(misses), proforma] + ) + "\n" + + # "w": the (single) writer process truncates on its first open of the + # run and appends afterwards, so a rerun cannot mix in stale results + queue.put((self.output_file, entries, False, "w")) + + def _build_proforma(self, prot_graph, inner_nodes, edges, peptide): + """ ProForma string for one peptide path. """ + # peptide index of the first residue of every sequence-bearing node + first_residue_index = {} + offset = 0 + for ni in inner_nodes: + first_residue_index[ni] = offset + offset += len(prot_graph.vs[ni]["aminoacid"]) + + n_term, c_term, by_residue = [], [], {} + for edge_id in edges: + edge = prot_graph.es[edge_id] + for key, delta in self._edge_mods(edge): + tag = self.mod_names.get(float(delta)) or self._delta_tag(delta) + if key in ("NPEPTERM", "NPROTERM"): + n_term.append(tag) + elif key in ("CPEPTERM", "CPROTERM"): + c_term.append(tag) + elif edge.target in first_residue_index: + by_residue.setdefault(first_residue_index[edge.target], []).append(tag) + # else: the edge points at a helper/terminal node and the key is + # not terminal — nothing to place (does not occur in generated + # graphs); never guess a position + + parts = [] + if n_term: + parts.append("".join("[{}]".format(t) for t in n_term) + "-") + for idx, aa in enumerate(peptide): + parts.append(aa) + parts.extend("[{}]".format(t) for t in by_residue.get(idx, ())) + if c_term: + parts.append("-" + "".join("[{}]".format(t) for t in c_term)) + return "".join(parts) + + def _edge_mods(self, edge): + """ Deduplicated (key, delta) modifications of one traversed edge. + + Or-wrapped qualifiers (collapsed parallel edges) are alternatives, not + a conjunction: only modifications present in EVERY branch are certain + for this traversal and get exported; branch-specific ones are dropped + with a one-time warning rather than fabricated onto every path. + """ + if "qualifiers" not in edge.attributes(): + return [] + mods = [] + for f in edge["qualifiers"] or []: + if isinstance(f, Or): + branches = [self._collect_mods(branch) for branch in f] + common = set.intersection(*map(set, branches)) if branches else set() + if any(set(b) - common for b in branches) and not self.warned_or_once: + self.warned_or_once = True + print( + "Warning: modifications differ between collapsed edge alternatives; " + "only modifications shared by all alternatives are exported. " + "Re-run with --no_collapsing_edges for exact per-path modifications." + ) + mods.extend(m for m in dict.fromkeys(x for b in branches for x in b) if m in common) + elif getattr(f, "type", None) in ("FIXMOD", "VARMOD"): + mods.append(self._feature_mod(f)) + return list(dict.fromkeys(mods)) + + def _collect_mods(self, qualifier): + out = [] + for f in qualifier or []: + if isinstance(f, Or): + for branch in f: + out.extend(self._collect_mods(branch)) + elif getattr(f, "type", None) in ("FIXMOD", "VARMOD"): + out.append(self._feature_mod(f)) + return out + + def _feature_mod(self, feature): + key, _, delta = feature.qualifiers["note"].rpartition(":") + return key, delta + + def _delta_tag(self, delta): + """ '79.966' -> '+79.966': explicitly signed, plain decimal notation. + + ProForma's mass grammar has no exponent and no non-finite forms, so + '1e-05' is rewritten as '+0.00001' and nan/inf are refused. + """ + if not isfinite(float(delta)): + raise ValueError( + "Cannot express the delta mass '{}' as a ProForma mass tag".format(delta) + ) + if "e" in delta.lower(): + delta = format(Decimal(delta), "f") + return delta if delta.startswith("-") else "+" + delta diff --git a/protgraph/protgraph.py b/protgraph/protgraph.py index 6547f15..a1db883 100644 --- a/protgraph/protgraph.py +++ b/protgraph/protgraph.py @@ -184,6 +184,7 @@ def __call__(self, parser, namespace, values, option_string=None): ("citus_peptide_export", cli.add_citus_peptide_export), ("sqlite_peptide_export", cli.add_sqlite_peptide_export), ("fasta_peptide_export", cli.add_fasta_peptide_export), + ("proforma_peptide_export", cli.add_proforma_peptide_export), ("trie_peptide_export", cli.add_trie_peptide_export), ("gremlin_graph_export", cli.add_gremlin_graph_export), ] diff --git a/tests/fixtures/minimal.txt b/tests/fixtures/minimal.txt new file mode 100644 index 0000000..d224351 --- /dev/null +++ b/tests/fixtures/minimal.txt @@ -0,0 +1,16 @@ +ID TEST_HUMAN Reviewed; 8 AA. +AC X9TEST; +DT 01-JAN-2020, integrated into UniProtKB/Swiss-Prot. +DT 01-JAN-2020, sequence version 1. +DT 01-JAN-2020, entry version 1. +DE RecName: Full=Minimal test protein; +GN Name=TEST; +OS Homo sapiens (Human). +OC Eukaryota; Metazoa; Chordata. +OX NCBI_TaxID=9606; +FT VARIANT 4 +FT /note="T -> M (test variant)" +FT /id="VAR_900001" +SQ SEQUENCE 8 AA; 920 MW; 0123456789ABCDEF CRC64; + MKCTMSAK +// diff --git a/tests/functional_test.py b/tests/functional_test.py index d27829c..6cfa6cf 100644 --- a/tests/functional_test.py +++ b/tests/functional_test.py @@ -214,3 +214,53 @@ def test_issue13(self): def test_issue41(self): args = protgraph.parse_args(["-n", "1", "-epepfasta", os.path.join(self.examples_path, "P49782.txt")]) protgraph.prot_graph(**args) + + def _run_proforma(self, extra_args, out): + """ ProForma export of tests/fixtures/minimal.txt (MKCTMSAK, VARIANT T4M). """ + fixture = os.path.abspath(os.path.join(os.path.dirname(__file__), "fixtures", "minimal.txt")) + args = protgraph.parse_args( + ["-epepproforma", "--pep_proforma_out", out, "--pep_hops", "10", "-d", "trypsin"] + + extra_args + self.procs_num + [fixture] + ) + protgraph.prot_graph(**args) + with open(out) as f: + return [ln.rstrip("\n").split("\t") for ln in f if ln.strip()] + + def test_export_pep_proforma_smoke(self): + args = protgraph.parse_args(["-epepproforma", "--pep_hops", "2"] + self.procs_num + self.example_files) + protgraph.prot_graph(**args) + + def test_export_pep_proforma_pins_mod_placement(self): + import tempfile + with tempfile.TemporaryDirectory() as tmp: + lines = self._run_proforma( + ["-ft", "VARIANT", "-fm", "C:57.021464", "-vm", "M:15.994915", + "--pep_proforma_mod_names", "15.994915=UNIMOD:35"], + os.path.join(tmp, "p.tsv"), + ) + assert all(len(x) == 5 for x in lines) + rows = {tuple(x) for x in lines} + proformas = {x[4] for x in lines} + # exact placement with coordinates: the fixed mod sits on the C + assert ("X9TEST", "3", "8", "0", "C[+57.021464]TMSAK") in rows + # a variable mod on a VARIANT-introduced residue (no reference + # position!) is placed on that residue, not dropped + assert "C[+57.021464]M[UNIMOD:35]MSAK" in proformas + # variable mods also export their unmodified counterparts + assert {"MK", "M[UNIMOD:35]K"} <= proformas + # the fixed mod is applied to every C, always signed + assert all("C[+57.021464]" in p for p in proformas if "C" in p) + + def test_export_pep_proforma_terminal_mods_and_overwrite(self): + import tempfile + with tempfile.TemporaryDirectory() as tmp: + out = os.path.join(tmp, "p.tsv") + lines = self._run_proforma(["-fm", "NPEPTERM:42.010565", "-vm", "CPEPTERM:79.966"], out) + proformas = {x[4] for x in lines} + assert all(p.startswith("[+42.010565]-") for p in proformas) + # variable C-terminal: modified and unmodified paths both exported + assert any(p.endswith("-[+79.966]") for p in proformas) + assert any(not p.endswith("-[+79.966]") for p in proformas) + # a rerun overwrites instead of appending stale results + again = self._run_proforma(["-fm", "NPEPTERM:42.010565", "-vm", "CPEPTERM:79.966"], out) + assert len(again) == len(lines) From 6610ad5c4c856b34584ff2f6e82e30ae963e5426 Mon Sep 17 00:00:00 2001 From: Dominik Lux Date: Fri, 4 Sep 2026 21:28:33 +0200 Subject: [PATCH 2/4] Apply review patch: qualifiers column, header row, merged-node fix, 0.3.13 Patch provided by the maintainer in the PR review: column order accession/proforma/start/end/misscleavages with an optional qualifiers column (--pep_proforma_write_qualifiers), a write-once TSV header, node-wise proforma assembly fixing modification placement inside merged nodes (EIN with modified I now yields EI[+x]N), README section moved after the FASTA walkthrough, version 0.3.13. --- README.md | 45 +++++++------- protgraph/cli.py | 8 ++- protgraph/export/peptides/pep_proforma.py | 74 ++++++++++++++--------- setup.py | 2 +- tests/functional_test.py | 7 ++- 5 files changed, 79 insertions(+), 57 deletions(-) diff --git a/README.md b/README.md index 01d1e31..c0ac51f 100644 --- a/README.md +++ b/README.md @@ -366,28 +366,6 @@ $ cat e_coli_with_selected_features_limited.fasta| grep "^>" | wc -l 648551 ``` -### Exporting peptides as ProForma 2.0 - -Peptides can also be exported in [ProForma 2.0](https://github.com/HUPO-PSI/ProForma) -notation (the HUPO-PSI standard for writing peptidoforms/proteoforms) into a -tab-separated file with the columns accession, start, end, misscleavages and -proforma: - -```shell -$ protgraph -d trypsin --pep_hops 3 -fm "C:57.021464" -vm "M:15.994915" -epepproforma --pep_proforma_out peptides.proforma.tsv examples/e_coli.dat -$ grep -m 2 "\[" peptides.proforma.tsv -P0A884 49 53 1 RC[+57.021464]HLR -P0A884 1 2 0 M[+15.994915]K -``` - -Modifications applied via `-fm`/`-vm` are known to ProtGraph only by their -delta mass, so they are written as mass-delta tags (`C[+57.021464]`), which is -valid ProForma. If you want named tags instead, declare a mapping for a delta -you passed: `--pep_proforma_mod_names "15.994915=UNIMOD:35"` then writes -`M[UNIMOD:35]`. Terminal modifications (`NPEPTERM` etc.) use the ProForma -terminal syntax (`[+42.010565]-PEPTIDE`). - - We can see that the number peptides within this FASTA is significantly reduced (to ~650 000 entries). **NOTE:** for setting an upper length-limit of peptides we need to set `-nm`. In case of not setting this parameter, longer peptides may be exported. Finally, since the FASTA-file still contains for some entries same sequences, we can concatenate these by using a different (and more sophisticated) FASTA-exporter: @@ -417,5 +395,26 @@ $ cat e_coli_compact.fasta | grep "^>" | wc -l Instead of directly generating a FASTA-file we first create a database, summarizing same sequences and headers. As a post-processing step, the database-entries are exported into FASTA. The generated FASTA has unique sequence as entries and offers some additionally insights. From the first entries we see peptides shared by exactly 2 proteins. Looking at the difference between the compact and non-compact FASTA, we see that ~244 000 entries could be summarized into already included entries in the FASTA. This generated FASTA-file can be used for identification. -**NOTE:** Protein-Graphs can contain large amounts of peptides/proteins. Do a dry run with the flags `-cnp`, `-cnpm` or `cnph` (or all of them) WITHOUT the export functionality first and examine the statistics output if it is feasible to generate a FASTA-file. Without a dry run it may happen that a protein like P04637 (P53 Human) with all possible peptides and variants is exported, which will +### Exporting peptides as ProForma 2.0 + +Peptides can also be exported in [ProForma 2.0](https://github.com/HUPO-PSI/ProForma) +notation (the HUPO-PSI standard for writing peptidoforms/proteoforms) into a +tab-separated file with the columns accession, start, end, misscleavages and +proforma: + +```shell +$ protgraph -d trypsin --pep_hops 3 -fm "C:57.021464" -vm "M:15.994915" -epepproforma --pep_proforma_out peptides.proforma.tsv examples/e_coli.dat +$ grep -m 2 "\[" peptides.proforma.tsv +P0A884 49 53 1 RC[+57.021464]HLR +P0A884 1 2 0 M[+15.994915]K +``` + +Modifications applied via `-fm`/`-vm` are known to ProtGraph only by their +delta mass, so they are written as mass-delta tags (`C[+57.021464]`), which is +valid ProForma. If you want named tags instead, declare a mapping for a delta +you passed: `--pep_proforma_mod_names "15.994915=UNIMOD:35"` then writes +`M[UNIMOD:35]`. Terminal modifications (`NPEPTERM` etc.) use the ProForma +terminal syntax (`[+42.010565]-PEPTIDE`). + +**NOTE:** Protein-Graphs can contain large amounts of peptides/proteins. Do a dry run with the flags `-cnp`, `-cnpm` or `cnph` (or all of them) WITHOUT the export functionality first and examine the statistics output if it is feasible to generate a peptide output file (e.g., fasta or proforma). Without a dry run it may happen that a protein like P04637 (P53 Human) with all possible peptides and variants is exported, which will very likely take up all your disk space. diff --git a/protgraph/cli.py b/protgraph/cli.py index 7c024cb..bbf2bbb 100644 --- a/protgraph/cli.py +++ b/protgraph/cli.py @@ -666,7 +666,7 @@ def add_proforma_peptide_export(group): group.add_argument( "--export_peptide_proforma", "-epepproforma", default=False, action="store_true", help="Set this flag to export peptides as ProForma 2.0 strings (PSI standard) " - "into a single tab-separated file (accession, start, end, misscleavages, proforma). " + "into a single tab-separated file (accession, proforma, start, end, misscleavages, qualifiers). " "Modifications from -fm/-vm are written as mass-delta tags (e.g. 'S[+79.966]')." ) group.add_argument( @@ -706,6 +706,12 @@ def _parse_mod_name(input: str): "of 'S[+79.966]'. Can be provided multiple times. No mass lookup is performed; " "the name is taken as declared." ) + group.add_argument( + "--pep_proforma_write_qualifiers", + action="store_true", + help="Adds the column 'qualifiers' into the ProForma table for tracing back how the peptide " + "was retrieved from the protein entry (e.g., 'VARIANT[22:22, T->L, REFERENCE]))." + ) def add_trie_peptide_export(group): diff --git a/protgraph/export/peptides/pep_proforma.py b/protgraph/export/peptides/pep_proforma.py index 588d421..ca14015 100644 --- a/protgraph/export/peptides/pep_proforma.py +++ b/protgraph/export/peptides/pep_proforma.py @@ -10,16 +10,15 @@ class PepProForma(PepFasta): Peptide exporter emitting ProForma 2.0 strings (HUPO-PSI standard notation for proteoforms and peptidoforms, LeDuc et al., J. Proteome Res. 2022). - Writes one tab-separated line per exported peptide path: - - accession start end misscleavages proforma + Writes one tab-separated line per exported peptide path with columns: + accession, proforma, start, end, misscleavages, qualifiers ProtGraph knows modifications only by their delta mass (-fm/-vm), so modifications are emitted as mass-delta tags by default ("S[+79.966]"), which is valid ProForma. A user-declared mapping (--pep_proforma_mod_names "79.966=UNIMOD:21") replaces a delta by a named tag ("S[UNIMOD:21]"). Terminal modifications (NPEPTERM etc.) use ProForma terminal syntax - ("[+42.010565]-PEPTIDE"). + ("[+42.010565]-PEPTIDE" or "[UNIMOD:1]-PEPTIDE"). A modification is placed on the residue its edge points at: FIXMOD/VARMOD features are attached to the in-edges of the (cloned) modified node when @@ -37,34 +36,37 @@ def start_up(self, **kwargs): self.output_file = kwargs["pep_proforma_out"] self.mod_names = dict(kwargs["pep_proforma_mod_names"] or []) self.warned_or_once = False + self.write_qualifiers = kwargs["pep_proforma_write_qualifiers"] + + def export(self, prot_graph, queue): + tsv_header = [ + "accession", "proforma", "start", "end", "misscleavages" + ] + (["qualifiers"] if self.write_qualifiers else []) + + queue.put((self.output_file, "\t".join(tsv_header) + "\n", True, "w")) + super().export(prot_graph, queue) def export_peptides(self, prot_graph, l_path_nodes, l_path_edges, l_peptide, l_miscleavages, queue): entries = "" - for peptide, nodes, edges, misses in zip(l_peptide, l_path_nodes, l_path_edges, l_miscleavages): + for nodes, edges, misses in zip(l_path_nodes, l_path_edges, l_miscleavages): # the first/last node carrying sequence; with terminal modifications # applied, nodes[1] can be an empty helper node (see annotate_ptms) - inner = [n for n in nodes[1:-1] if prot_graph.vs[n]["aminoacid"]] - acc = self._get_accession_or_isoform(prot_graph.vs[inner[0]]) - start_pos = self._get_position_or_isoform_position(prot_graph.vs[inner[0]]) - end_pos = self._get_position_or_isoform_position(prot_graph.vs[inner[-1]], end=True) - proforma = self._build_proforma(prot_graph, inner, edges, peptide) + acc = self._get_accession_or_isoform(prot_graph.vs[nodes[1]]) + start_pos = self._get_position_or_isoform_position(prot_graph.vs[nodes[1]]) + end_pos = self._get_position_or_isoform_position(prot_graph.vs[nodes[-2]], end=True) + proforma = self._build_proforma(prot_graph, nodes[1:-1], edges) + l_str_qualifiers = self._get_qualifiers(prot_graph, edges) entries += "\t".join( - [acc, str(start_pos), str(end_pos), str(misses), proforma] + [acc, proforma, str(start_pos), str(end_pos), str(misses)] + + (l_str_qualifiers if self.write_qualifiers else []) ) + "\n" - # "w": the (single) writer process truncates on its first open of the - # run and appends afterwards, so a rerun cannot mix in stale results + # w: ensures to have no stale entries, overwriting the file, if it exists. queue.put((self.output_file, entries, False, "w")) - def _build_proforma(self, prot_graph, inner_nodes, edges, peptide): + def _build_proforma(self, prot_graph, inner_nodes, edges): """ ProForma string for one peptide path. """ - # peptide index of the first residue of every sequence-bearing node - first_residue_index = {} - offset = 0 - for ni in inner_nodes: - first_residue_index[ni] = offset - offset += len(prot_graph.vs[ni]["aminoacid"]) - + # Get all modifications by iterating over the edges n_term, c_term, by_residue = [], [], {} for edge_id in edges: edge = prot_graph.es[edge_id] @@ -74,18 +76,32 @@ def _build_proforma(self, prot_graph, inner_nodes, edges, peptide): n_term.append(tag) elif key in ("CPEPTERM", "CPROTERM"): c_term.append(tag) - elif edge.target in first_residue_index: - by_residue.setdefault(first_residue_index[edge.target], []).append(tag) - # else: the edge points at a helper/terminal node and the key is - # not terminal — nothing to place (does not occur in generated - # graphs); never guess a position + else: + by_residue.setdefault(edge.target, []).append((key, tag)) + # Build the final proforma sequence parts = [] if n_term: parts.append("".join("[{}]".format(t) for t in n_term) + "-") - for idx, aa in enumerate(peptide): - parts.append(aa) - parts.extend("[{}]".format(t) for t in by_residue.get(idx, ())) + for ie in inner_nodes: + aas = prot_graph.vs[ie]["aminoacid"] + if len(aas) == 1: + # Simply append UniMod + parts.append(aas) + parts.extend("[{}]".format(t[1]) for t in by_residue.get(ie, ())) + else: + # We need to add PTM on correct position in merged node: + # E.G. Merged node "EIN", Modification I->+23.123 + # should yield "EI[+23.123]N" and NOT "EIN[+23.123]"! + by_aa = {} + for mods in by_residue.values(): + for aa, mod in mods: + by_aa.setdefault(aa, []).append(mod) + for aa, mod in by_aa.items(): + # We can savely use replace as this case only can happen in FIXMODs + aas = aas.replace(aa, aa + "".join("[{}]".format(t) for t in mod)) + parts.append(aas) + if c_term: parts.append("-" + "".join("[{}]".format(t) for t in c_term)) return "".join(parts) diff --git a/setup.py b/setup.py index 79ae5ed..49ab675 100644 --- a/setup.py +++ b/setup.py @@ -11,7 +11,7 @@ setup( name='protgraph', - version='0.3.12', + version='0.3.13', author="Dominik Lux", description="ProtGraph, a graph generator for proteins.", long_description=long_description, diff --git a/tests/functional_test.py b/tests/functional_test.py index 6cfa6cf..b09453b 100644 --- a/tests/functional_test.py +++ b/tests/functional_test.py @@ -240,9 +240,9 @@ def test_export_pep_proforma_pins_mod_placement(self): ) assert all(len(x) == 5 for x in lines) rows = {tuple(x) for x in lines} - proformas = {x[4] for x in lines} + proformas = {x[1] for x in lines} # exact placement with coordinates: the fixed mod sits on the C - assert ("X9TEST", "3", "8", "0", "C[+57.021464]TMSAK") in rows + assert ("X9TEST", "C[+57.021464]TMSAK", "3", "8", "0",) in rows # a variable mod on a VARIANT-introduced residue (no reference # position!) is placed on that residue, not dropped assert "C[+57.021464]M[UNIMOD:35]MSAK" in proformas @@ -256,7 +256,8 @@ def test_export_pep_proforma_terminal_mods_and_overwrite(self): with tempfile.TemporaryDirectory() as tmp: out = os.path.join(tmp, "p.tsv") lines = self._run_proforma(["-fm", "NPEPTERM:42.010565", "-vm", "CPEPTERM:79.966"], out) - proformas = {x[4] for x in lines} + proformas = {x[1] for x in lines} + proformas.discard("proforma") assert all(p.startswith("[+42.010565]-") for p in proformas) # variable C-terminal: modified and unmodified paths both exported assert any(p.endswith("-[+79.966]") for p in proformas) From c8538502335a6f39d10fb91f7ac6b5902c2ced29 Mon Sep 17 00:00:00 2001 From: Oliver Kohlbacher Date: Fri, 4 Sep 2026 21:28:33 +0200 Subject: [PATCH 3/4] Take start/end from sequence-bearing nodes With a terminal modification applied, nodes[1] / nodes[-2] can be the empty helper node whose position 0 does not exist in 1-based protein coordinates; an NPEPTERM peptide reported start=0. Positions now come from the first and last node carrying sequence. Co-Authored-By: Claude Fable 5 --- protgraph/export/peptides/pep_proforma.py | 14 +++++++++----- 1 file changed, 9 insertions(+), 5 deletions(-) diff --git a/protgraph/export/peptides/pep_proforma.py b/protgraph/export/peptides/pep_proforma.py index ca14015..3a885fb 100644 --- a/protgraph/export/peptides/pep_proforma.py +++ b/protgraph/export/peptides/pep_proforma.py @@ -49,11 +49,15 @@ def export(self, prot_graph, queue): def export_peptides(self, prot_graph, l_path_nodes, l_path_edges, l_peptide, l_miscleavages, queue): entries = "" for nodes, edges, misses in zip(l_path_nodes, l_path_edges, l_miscleavages): - # the first/last node carrying sequence; with terminal modifications - # applied, nodes[1] can be an empty helper node (see annotate_ptms) - acc = self._get_accession_or_isoform(prot_graph.vs[nodes[1]]) - start_pos = self._get_position_or_isoform_position(prot_graph.vs[nodes[1]]) - end_pos = self._get_position_or_isoform_position(prot_graph.vs[nodes[-2]], end=True) + # start/end come from the first/last node carrying sequence: with + # terminal modifications applied, nodes[1] (or nodes[-2]) can be an + # empty helper node (see annotate_ptms) whose position 0 does not + # exist in 1-based protein coordinates + first = next(n for n in nodes[1:-1] if prot_graph.vs[n]["aminoacid"]) + last = next(n for n in nodes[-2:0:-1] if prot_graph.vs[n]["aminoacid"]) + acc = self._get_accession_or_isoform(prot_graph.vs[first]) + start_pos = self._get_position_or_isoform_position(prot_graph.vs[first]) + end_pos = self._get_position_or_isoform_position(prot_graph.vs[last], end=True) proforma = self._build_proforma(prot_graph, nodes[1:-1], edges) l_str_qualifiers = self._get_qualifiers(prot_graph, edges) entries += "\t".join( From 09a6a3fb5e34cbe631c18e9bc5acd002c899d70a Mon Sep 17 00:00:00 2001 From: Oliver Kohlbacher Date: Fri, 4 Sep 2026 23:32:49 +0200 Subject: [PATCH 4/4] Place merged-node mods per node and position; qualifiers as one column MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Full-scale end-to-end testing surfaced two defects in the previous state: Merged-node placement pooled (letter, tag) pairs across ALL nodes of the path and applied them with sequential str.replace, so a peptide spanning several nodes that each carry the same-letter modification received the accumulated tags on every matching residue (silently inflating the mass, e.g. C[+57][+57] under fixed carbamidomethyl with one miscleavage), and a residue letter occurring inside an already-written tag's text got tags injected mid-tag (C[UNIM[+15.99]OD:4]). Placement now uses only the node's own modifications, locating each one via the feature's position relative to the node start (falling back to the first matching residue for features without a usable position), and assembles the string per residue index — no replace, each modification exactly once. The qualifiers option wrote one tab-separated FIELD per qualifier, giving a ragged TSV (5-13 columns under a 6-column header); the qualifier strings are now comma-joined into the single declared column, present (possibly empty) on every row. New fixture tests/fixtures/merged.txt pins both: EI[+23.123]NK for an interior modification of a merged node, single tags on a two-node miscleaved peptide, no adjacent duplicate and no nested tags, and a uniform 6-column table with the flag set. Co-Authored-By: Claude Fable 5 --- protgraph/export/peptides/pep_proforma.py | 61 ++++++++++++++++------- tests/fixtures/merged.txt | 12 +++++ tests/functional_test.py | 47 +++++++++++++++++ 3 files changed, 101 insertions(+), 19 deletions(-) create mode 100644 tests/fixtures/merged.txt diff --git a/protgraph/export/peptides/pep_proforma.py b/protgraph/export/peptides/pep_proforma.py index 3a885fb..1196ea8 100644 --- a/protgraph/export/peptides/pep_proforma.py +++ b/protgraph/export/peptides/pep_proforma.py @@ -1,6 +1,8 @@ from decimal import Decimal from math import isfinite +from Bio.SeqFeature import UnknownPosition + from protgraph.export.peptides.pep_fasta import PepFasta from protgraph.graph_collapse_edges import Or @@ -59,10 +61,13 @@ def export_peptides(self, prot_graph, l_path_nodes, l_path_edges, l_peptide, l_m start_pos = self._get_position_or_isoform_position(prot_graph.vs[first]) end_pos = self._get_position_or_isoform_position(prot_graph.vs[last], end=True) proforma = self._build_proforma(prot_graph, nodes[1:-1], edges) - l_str_qualifiers = self._get_qualifiers(prot_graph, edges) entries += "\t".join( [acc, proforma, str(start_pos), str(end_pos), str(misses)] - + (l_str_qualifiers if self.write_qualifiers else []) + # one 'qualifiers' COLUMN, matching the header: the individual + # qualifier strings are comma-joined (as in the FASTA headers), + # and the field is present, possibly empty, on every row + + ([",".join(self._get_qualifiers(prot_graph, edges))] + if self.write_qualifiers else []) ) + "\n" # w: ensures to have no stale entries, overwriting the file, if it exists. @@ -74,14 +79,14 @@ def _build_proforma(self, prot_graph, inner_nodes, edges): n_term, c_term, by_residue = [], [], {} for edge_id in edges: edge = prot_graph.es[edge_id] - for key, delta in self._edge_mods(edge): + for key, delta, loc_end in self._edge_mods(edge): tag = self.mod_names.get(float(delta)) or self._delta_tag(delta) if key in ("NPEPTERM", "NPROTERM"): n_term.append(tag) elif key in ("CPEPTERM", "CPROTERM"): c_term.append(tag) else: - by_residue.setdefault(edge.target, []).append((key, tag)) + by_residue.setdefault(edge.target, []).append((key, tag, loc_end)) # Build the final proforma sequence parts = [] @@ -89,29 +94,46 @@ def _build_proforma(self, prot_graph, inner_nodes, edges): parts.append("".join("[{}]".format(t) for t in n_term) + "-") for ie in inner_nodes: aas = prot_graph.vs[ie]["aminoacid"] - if len(aas) == 1: + node_mods = by_residue.get(ie, ()) + if not node_mods: + parts.append(aas) + elif len(aas) == 1: # Simply append UniMod parts.append(aas) - parts.extend("[{}]".format(t[1]) for t in by_residue.get(ie, ())) + parts.extend("[{}]".format(t) for _, t, _ in node_mods) else: - # We need to add PTM on correct position in merged node: - # E.G. Merged node "EIN", Modification I->+23.123 - # should yield "EI[+23.123]N" and NOT "EIN[+23.123]"! - by_aa = {} - for mods in by_residue.values(): - for aa, mod in mods: - by_aa.setdefault(aa, []).append(mod) - for aa, mod in by_aa.items(): - # We can savely use replace as this case only can happen in FIXMODs - aas = aas.replace(aa, aa + "".join("[{}]".format(t) for t in mod)) - parts.append(aas) + # We need to add each PTM on its correct position in the merged + # node: e.g. merged node "EIN", modification I->+23.123 should + # yield "EI[+23.123]N" and NOT "EIN[+23.123]". Only THIS node's + # modifications are placed (pooling across the whole path + # duplicated same-letter tags), each exactly once, located via + # the feature's position relative to the node start; a feature + # without a usable position (e.g. on a variant residue) falls + # back to the first matching residue letter. + node_pos = self._get_position_or_isoform_position(prot_graph.vs[ie]) + tags_at = {} + for key, tag, loc_end in node_mods: + idx = None + if loc_end is not None and isinstance(node_pos, int): + offset = loc_end - node_pos + if 0 <= offset < len(aas): + idx = offset + if idx is None: + idx = aas.find(key[-1]) # positional keys end in the residue letter + if idx == -1: + idx = 0 # never drop a modification the path carries + tags_at.setdefault(idx, []).append(tag) + parts.append("".join( + aa + "".join("[{}]".format(t) for t in tags_at.get(i, ())) + for i, aa in enumerate(aas) + )) if c_term: parts.append("-" + "".join("[{}]".format(t) for t in c_term)) return "".join(parts) def _edge_mods(self, edge): - """ Deduplicated (key, delta) modifications of one traversed edge. + """ Deduplicated (key, delta, position) modifications of one traversed edge. Or-wrapped qualifiers (collapsed parallel edges) are alternatives, not a conjunction: only modifications present in EVERY branch are certain @@ -149,7 +171,8 @@ def _collect_mods(self, qualifier): def _feature_mod(self, feature): key, _, delta = feature.qualifiers["note"].rpartition(":") - return key, delta + end = feature.location.end + return key, delta, None if isinstance(end, UnknownPosition) else int(end) def _delta_tag(self, delta): """ '79.966' -> '+79.966': explicitly signed, plain decimal notation. diff --git a/tests/fixtures/merged.txt b/tests/fixtures/merged.txt new file mode 100644 index 0000000..0da3446 --- /dev/null +++ b/tests/fixtures/merged.txt @@ -0,0 +1,12 @@ +ID MERGED_TEST Reviewed; 11 AA. +AC X9MRG0; +DT 01-JAN-2020, integrated into UniProtKB/Swiss-Prot. +DT 01-JAN-2020, sequence version 1. +DT 01-JAN-2020, entry version 1. +DE RecName: Full=Merged node modification test; +OS Homo sapiens (Human). +OC Eukaryota; Metazoa; Chordata. +OX NCBI_TaxID=9606; +SQ SEQUENCE 11 AA; 1200 MW; 0123456789ABCDEF CRC64; + KEINKACKACK +// diff --git a/tests/functional_test.py b/tests/functional_test.py index b09453b..d80c3fe 100644 --- a/tests/functional_test.py +++ b/tests/functional_test.py @@ -265,3 +265,50 @@ def test_export_pep_proforma_terminal_mods_and_overwrite(self): # a rerun overwrites instead of appending stale results again = self._run_proforma(["-fm", "NPEPTERM:42.010565", "-vm", "CPEPTERM:79.966"], out) assert len(again) == len(lines) + + def test_export_pep_proforma_places_mods_inside_merged_nodes_exactly_once(self): + import re + import tempfile + fixture = os.path.abspath(os.path.join(os.path.dirname(__file__), "fixtures", "merged.txt")) + with tempfile.TemporaryDirectory() as tmp: + out = os.path.join(tmp, "p.tsv") + args = protgraph.parse_args( + [ + "-epepproforma", "--pep_proforma_out", out, "--pep_hops", "10", + "-d", "trypsin", "--pep_miscleavages", "2", + "-fm", "I:23.123", "-fm", "C:57.021464", + ] + self.procs_num + [fixture] + ) + protgraph.prot_graph(**args) + with open(out) as f: + next(f) # header + proformas = {ln.rstrip("\n").split("\t")[1] for ln in f if ln.strip()} + # a mod on an interior residue of a merged node sits on that residue + assert "EI[+23.123]NK" in proformas + # a miscleaved peptide spanning two nodes that EACH carry the same + # fixed mod gets one tag per residue, not the pooled duplicates + assert "AC[+57.021464]KAC[+57.021464]K" in proformas + for p in proformas: + assert "[+57.021464][+57.021464]" not in p and "[+23.123][+23.123]" not in p + assert not re.search(r"\[[^]]*\[", p), p # no tag nested inside a tag + + def test_export_pep_proforma_qualifiers_is_one_column(self): + import tempfile + fixture = os.path.abspath(os.path.join(os.path.dirname(__file__), "fixtures", "merged.txt")) + with tempfile.TemporaryDirectory() as tmp: + out = os.path.join(tmp, "p.tsv") + args = protgraph.parse_args( + [ + "-epepproforma", "--pep_proforma_out", out, "--pep_hops", "10", + "-d", "trypsin", "--pep_miscleavages", "2", + "-fm", "C:57.021464", "--pep_proforma_write_qualifiers", + ] + self.procs_num + [fixture] + ) + protgraph.prot_graph(**args) + with open(out) as f: + widths = {len(ln.rstrip("\n").split("\t")) for ln in f if ln.strip()} + f.seek(0) + next(f) + quals = [ln.rstrip("\n").split("\t")[5] for ln in f if ln.strip()] + assert widths == {6} # header and every row: exactly one qualifiers column + assert any("FIXMOD" in q for q in quals)