Metadata-Version: 2.5
Name: genotopep
Version: 0.2.0
Summary: Reproducible peptide candidate generation from assemblies and predicted proteins
Project-URL: Homepage, https://github.com/ljunwon1114/GenoToPep
Project-URL: Repository, https://github.com/ljunwon1114/GenoToPep
Project-URL: Issues, https://github.com/ljunwon1114/GenoToPep/issues
Author: Jun Won Lee
License-Expression: GPL-3.0-or-later
License-File: LICENSE
Keywords: bioinformatics,genomics,metagenomics,peptides
Classifier: Development Status :: 3 - Alpha
Classifier: Environment :: Console
Classifier: Intended Audience :: Science/Research
Classifier: License :: OSI Approved :: GNU General Public License v3 or later (GPLv3+)
Classifier: Operating System :: OS Independent
Classifier: Programming Language :: Python :: 3
Classifier: Programming Language :: Python :: 3.10
Classifier: Programming Language :: Python :: 3.11
Classifier: Programming Language :: Python :: 3.12
Classifier: Topic :: Scientific/Engineering :: Bio-Informatics
Requires-Python: >=3.10
Requires-Dist: pyrodigal<4,>=3.6
Requires-Dist: pyyaml<7,>=6
Provides-Extra: dev
Requires-Dist: build>=1.2; extra == 'dev'
Requires-Dist: pytest<10,>=8; extra == 'dev'
Requires-Dist: ruff<1,>=0.9; extra == 'dev'
Description-Content-Type: text/markdown

# GenoToPep

GenoToPep is a reproducible command-line package for generating peptide candidate universes from metagenome assemblies or predicted-protein FASTA files.

## Installation

Python 3.10 or newer is required.

### PyPI

```bash
python -m pip install genotopep
genotopep --version
genotopep --help
```

### Bioconda

GenoToPep is available from Bioconda. The current public Bioconda release is `0.1.0.1` (`noarch`, build `pyhdfd78af_0`):

```bash
conda install -c conda-forge -c bioconda genotopep
```

For an exact reproducible install of the current public release:

```bash
conda install -c conda-forge -c bioconda genotopep=0.1.0.1
```

The local 0.2.0 development tree documented below is not yet published to Bioconda.

## Input routes

- Nucleotide assembly: `--assembly sample.fna[.gz]`
- Predicted proteins: `--proteins sample.faa[.gz]`
- Predicted proteins with genomic provenance: `--proteins sample.faa --gff sample.gff3[.gz]`

`--assembly` and `--proteins` are mutually exclusive. GFF3 is optional and is only valid with FAA input.

### FNA ORF calling

Assembly input is translated with Pyrodigal metagenome mode [1,2]. The default short-ORF threshold is:

```text
min_protein_aa = 10
min_gene_nt = 10 * 3 + 3 = 33 nt
closed = false
max_overlap = 0
```

Pyrodigal still applies its coding and start-site models: `33 nt` is an eligible minimum, not a guarantee that every open reading frame of that length is emitted. Contig-edge partial calls are retained and marked in the ORF provenance output. Terminal stop symbols are not written to protein FASTA.

## Generation strategies

### Intact

Emits the complete supplied or Pyrodigal-called protein when its length is within the inclusive range.

```bash
genotopep generate intact \
  --assembly assembly.fna.gz \
  --min-length 10 --max-length 100 \
  --output results/assembly-intact
```

```bash
genotopep generate intact \
  --proteins proteins.faa.gz \
  --min-length 10 --max-length 100 \
  --output results/proteins-intact
```

### Windows

Emits contiguous subsequences for every inclusive length, using stride 1 by default.

```bash
genotopep generate windows \
  --assembly assembly.fna.gz \
  --min-length 10 --max-length 50 --stride 1 \
  --output results/assembly-windows
```

```bash
genotopep generate windows \
  --proteins proteins.faa.gz \
  --min-length 10 --max-length 50 --stride 1 \
  --output results/proteins-windows
```

### Cleavage

Builds a boundary union from selected deterministic rules, creates elementary fragments, and emits every length-valid joining of adjacent fragments.

```bash
genotopep generate cleavage \
  --assembly assembly.fna.gz \
  --cleavage-mode all \
  --min-length 10 --max-length 100 \
  --output results/assembly-cleavage
```

```bash
genotopep generate cleavage \
  --proteins proteins.faa.gz \
  --cleavage-mode trypsin cnbr \
  --min-length 10 --max-length 100 \
  --output results/proteins-cleavage
```

Supported lowercase modes follow simplified deterministic profiles derived from published specificity summaries and ExPASy PeptideCutter conventions [3-5]:

- `trypsin`: simplified deterministic profile that cuts after K/R unless followed by P; reported context-dependent exceptions and counter-exceptions are not modeled [3,5]
- `chymotrypsin`: high-specificity profile that cuts after F/Y unless followed by P and after W unless followed by M/P [3,5]
- `cnbr`: boundary after M; product chemistry, methionine oxidation, and condition-dependent incomplete cleavage are not modeled [5,6]
- `all`: union of trypsin, chymotrypsin, and CNBr boundaries

Omitting `--cleavage-mode` resolves to `all`. Multiple individual modes are accepted. `all` cannot be combined with an individual mode. The cited sources support residue-specific boundary conventions; the union of selected boundaries and enumeration of every contiguous adjacent-fragment joining are GenoToPep-defined candidate-generation semantics.

## Pre-generation estimates

`estimate` reports the exact occurrence count needed for workload planning without writing peptide outputs or performing sequence deduplication. An occurrence is a location-derived candidate, including duplicate sequences: if the same sequence is generated from multiple proteins or coordinates, each location is counted separately. Therefore `unique peptides <= retained occurrences`; the retained occurrence count is an exact upper bound on the final exact-sequence unique peptide count. Runtime, occurrence-TSV rows, and candidate-processing work scale primarily with occurrences rather than unique sequences.

The exclusion count is occurrence-level too. Because the exclusion TSV aggregates by source protein, strategy, and reason, one exclusion row can represent many excluded occurrences. `counts.exclusions` is the excluded-candidate count, whereas `counts.exclusion_groups` is the number of written aggregate rows.

```bash
genotopep estimate windows \
  --proteins proteins.faa.gz \
  --min-length 10 --max-length 50 --stride 1
```

```bash
genotopep estimate intact cleavage \
  --assembly assembly.fna.gz \
  --cleavage-mode trypsin cnbr \
  --min-length 10 --max-length 100
```

```bash
genotopep estimate \
  --assembly assembly.fna.gz \
  --cleavage-mode all \
  --min-length 10 --max-length 100 --stride 1
```

Omitting strategy selection defaults to `all`; explicitly specifying strategy `all` has the same effect. Both expand to `intact`, `windows`, and `cleavage`. `all` cannot be combined with individual strategy names. For assembly input, Pyrodigal ORF calling is performed to obtain the protein sequences required for exact counting, but ORF and peptide output files are not written. Reported occurrences are exact pre-deduplication counts and an upper bound on unique peptides. The exact unique count is produced by generation and exact-sequence deduplication.

Estimate JSON schema version 1 includes strategy-specific counts, arithmetic totals across the selected strategies, and non-blocking scale warnings. A strategy with at least 10,000,000 candidate events (`occurrences + exclusions`) receives a `large_candidate_universe` warning. The threshold does not limit generation, require `--force`, or predict runtime, memory, compressed size, or temporary-database size. Exhaustive windows grow rapidly: one canonical 300-aa protein produces 11,111 windows for lengths 10-50 at stride 1, so one million such proteins would produce 11,111,000,000 occurrences. Run `estimate` before large generation jobs and provision compute and storage from representative local measurements.

## Changes from 0.1.x

GenoToPep 0.2.0 is an output-schema migration. Manifest schema 2 is not resume-compatible with schema-1 or 0.1.x output directories. Exclusion TSV files now aggregate rows by source and reason while preserving exact excluded-occurrence totals in `excluded_occurrences`; they no longer enumerate every excluded candidate sequence. Unique peptide FASTA records are ordered lexicographically by exact amino-acid sequence rather than first occurrence. Estimate JSON schema 1 is a separate contract from manifest schema 2. Re-run generation under 0.2.0 when schema-2 outputs are required; do not merge 0.1.x and 0.2.0 partial outputs.

## Outputs

For input `sample.proteins.faa.gz` and strategy `windows`, the main outputs are:

```text
sample.proteins.windows.peptides.faa.gz
sample.proteins.windows.peptide_occurrences.tsv.gz
sample.proteins.windows.exclusions.tsv.gz
sample.proteins.command.txt
sample.proteins.run-config.yaml
sample.proteins.run.log
sample.proteins.manifest.json
sample.proteins.checksums.sha256
```

FNA runs additionally write:

```text
sample.orfs.faa.gz
sample.orfs.fna.gz
sample.orfs.gff3.gz
sample.orfs.tsv.gz
```

The ORF TSV `protein_sequence` and `nucleotide_sequence` columns intentionally duplicate the sequences in `orfs.faa.gz` and `orfs.fna.gz`, respectively, so each ORF row is self-contained for tabular analysis. GenoToPep 0.2.0 does not provide an option to omit these columns; account for this duplication when planning storage for very large assemblies.

### Relationships among output files

- `*.peptides.faa.gz` contains one record per exact peptide sequence; the FASTA header is the unique `peptide_id` (`>sha256:...`).
- `*.peptide_occurrences.tsv.gz` joins to the peptide FASTA on `peptide_id`. Many occurrence rows can reference one FASTA record, so this is a many-to-one relationship from occurrences to peptides.
- For FNA input, `*.orfs.tsv.gz` joins `id` to `*.peptide_occurrences.tsv.gz` `source_gene_id`. FAA input does not produce an ORF TSV; `source_gene_id` then refers to the supplied protein identifier.
- `*.exclusions.tsv.gz` has no `peptide_id`. For FNA input its `source_gene_id` can join to ORF `id`, but exclusions and retained occurrence rows are mutually exclusive at the candidate-event level.

The occurrence TSV has one data row per retained candidate occurrence. The exclusion TSV instead aggregates excluded candidates by source protein, strategy, and reason: `counts.exclusions` equals the sum of `excluded_occurrences`, while `counts.exclusion_groups` equals the exclusion TSV data-row count. Occurrence and exclusion TSV rows must not be added as row counts. Use `counts.occurrences + counts.exclusions` for the complete candidate-event total; aggregated exclusions do not retain each excluded candidate sequence or coordinate.

The peptide FASTA contains one record per exact amino-acid sequence. Its ID is `sha256:<digest>`. Exact sequence, rather than the digest, is the internal deduplication key. The occurrence TSV retains every source, protein-relative coordinate, genomic coordinate when available, strategy, window, and cleavage occurrence.

### Occurrence TSV columns

Interpret `coordinate_scope` before the genomic-coordinate columns. `protein_and_genomic` means columns `parent_contig` through `peptide_end_nt_1based` are populated from an FNA-called ORF or a successful FAA+GFF3 mapping. `protein_only` means genomic mapping is unavailable; for example, an empty `gene_strand` in an FAA run without GFF3 means unknown genomic context, not a strandless feature.

The remaining conditional groups are:

- `window_length` and `window_stride`: populated only for `windows`; otherwise not applicable to that generation strategy.
- `cleavage_modes`, `elementary_fragment_count`, and `internal_boundaries_joined`: populated only for `cleavage`; otherwise not applicable to that generation strategy.
- `orf_caller`, `orf_caller_version`, and `min_gene_nt`: populated only for FNA input; not applicable to supplied FAA input.

All integer coordinates are inclusive and 1-based. Empty values are serialized as zero-length TSV fields, not the strings `None` or `null`. An empty genomic field with `coordinate_scope=protein_only` means unavailable because no genomic mapping was supplied; an empty strategy- or input-specific field means not applicable. Readers such as pandas normally convert both cases to `NaN`, so retain `coordinate_scope`, `generation_strategy`, and `source_input_type` when interpreting nulls.

| Column | Type | Definition and empty-value meaning |
|---|---|---|
| `peptide_id` | string | `sha256:<digest>` of the exact normalized peptide sequence |
| `peptide_sequence` | string | Exact candidate amino-acid sequence |
| `peptide_length` | integer | Candidate length in amino acids |
| `generation_strategy` | enum | `intact`, `windows`, or `cleavage` |
| `source_input_type` | enum | `faa` or `fna` |
| `source_gene_id` | string | Source protein or called-ORF identifier |
| `parent_protein_id` | string | Parent protein identifier; currently equal to `source_gene_id` |
| `original_protein_header` | string | Original FAA header, or called-ORF identifier for FNA input |
| `protein_length` | integer | Normalized parent-protein length in amino acids |
| `peptide_start_aa_1based` | integer | Peptide start in the parent protein |
| `peptide_end_aa_1based` | integer | Peptide end in the parent protein |
| `coordinate_scope` | enum | Interpret first: `protein_and_genomic` permits genomic-column use; `protein_only` means genomic mapping is unavailable |
| `parent_contig` | string or empty | Parent contig; empty for FAA input without GFF3 |
| `gene_start_nt_1based` | integer or empty | CDS/ORF genomic start; empty for FAA input without GFF3 |
| `gene_end_nt_1based` | integer or empty | CDS/ORF genomic end; empty for FAA input without GFF3 |
| `gene_strand` | enum or empty | `+` or `-`; empty for FAA input without GFF3 because strand is unknown |
| `peptide_start_nt_1based` | integer or empty | Peptide genomic start; empty when genomic mapping is unavailable |
| `peptide_end_nt_1based` | integer or empty | Peptide genomic end; empty when genomic mapping is unavailable |
| `window_length` | integer or empty | Window length for `windows`; empty otherwise |
| `window_stride` | integer or empty | Window stride for `windows`; empty otherwise |
| `cleavage_modes` | comma-separated string or empty | Resolved modes for `cleavage`; empty otherwise |
| `elementary_fragment_count` | integer or empty | Adjacent elementary fragments joined for `cleavage`; empty otherwise |
| `internal_boundaries_joined` | integer or empty | Internal cleavage boundaries crossed for `cleavage`; empty otherwise |
| `orf_caller` | string or empty | ORF caller for FNA input; empty for FAA input |
| `orf_caller_version` | string or empty | ORF-caller version for FNA input; empty for FAA input |
| `min_gene_nt` | integer or empty | Minimum gene length for FNA ORF calling; empty for FAA input |

### Exclusion TSV columns

Noncanonical exclusions are aggregated to one row per source protein, strategy, and reason. `counts.exclusions` remains the exact number of excluded candidate occurrences; `counts.exclusion_groups` is the number of written exclusion rows.

| Column | Definition |
|---|---|
| `source_input_type` | `faa` or `fna` |
| `source_gene_id` | Source protein or called-ORF identifier |
| `parent_protein_id` | Parent protein identifier |
| `original_protein_header` | Original FAA header or called-ORF identifier |
| `protein_length` | Normalized parent-protein length |
| `generation_strategy` | Strategy whose candidates were excluded |
| `reason` | Exclusion reason; currently `noncanonical_amino_acid` |
| `excluded_occurrences` | Exact number of excluded candidate occurrences in the group |
| `noncanonical_positions_aa_1based` | Comma-separated 1-based noncanonical residue positions in the parent protein |
| `noncanonical_residues` | Residues corresponding to the reported positions |

Only the canonical 20 amino acids are admitted to predictor-facing peptide FASTA. Candidates containing other symbols are written to the exclusions table with a reason.

### Input normalization and identifiers

FAA sequences are converted to uppercase and trailing terminal `*` symbols are removed before length filtering, generation, peptide hashing, and output. FNA sequences and Pyrodigal protein translations are also converted to uppercase. The run configuration records `normalization.case: upper`. Internal tab, control, and DEL characters in FASTA headers are rejected because they are not safely interoperable with tabular command-line tools. FASTA identifiers must be unique within each input file; duplicates are hard errors.

### Integrity, payload identity, and resume

`checksums.sha256` covers the complete run inventory, including timestamped provenance files. Manifest schema version 2 additionally records role-specific SHA-256 values over the exact compressed bytes of the three deterministic scientific payload files—peptide FASTA, occurrence TSV, and exclusion TSV—and a combined `payload_sha256`. Scientific payloads use gzip `compresslevel 6`, an empty embedded filename, and `mtime 0`; these settings are recorded in `run-config.yaml`. The combined digest processes roles in lexical order using `UTF-8 role + NUL + lowercase ASCII SHA-256 hex + LF` framing; it excludes output-directory paths, timestamps, and command text, so equivalent payloads can be compared across output directories under the same supported compression environment.

Compressed-byte hashes can depend on the Python/zlib compression implementation. For cross-environment scientific comparison, compare decompressed payload content in addition to the recorded exact compressed-byte hashes.

These self-authored hashes detect accidental corruption and verify internal consistency; they are not authentication against deliberate modification. Anyone who can rewrite the complete output directory can replace payloads and recompute every manifest and checksum. Adversarial tamper resistance requires an external trust anchor such as a separately stored expected digest, signed manifest, or immutable repository record.

Outputs are written through an incomplete sibling directory and renamed only after successful completion. Failed runs close the temporary SQLite database and remove their incomplete directory. Concurrent GenoToPep generation attempts for the same destination are serialized by an exclusive sibling lock file; existing outputs are checked again before publication and are not intentionally overwritten. A hard process kill can bypass Python cleanup and leave a stale lock or `.incomplete.<UUID>` directory. After confirming that no matching generation process is running, remove only the stale artifacts belonging to the terminated invocation. The lock coordinates GenoToPep invocations; it cannot prevent an unrelated external program from racing to create the same path. `--resume` accepts only an identical, complete run after schema, code identity, input type, input/GFF SHA-256, parameters, normalization, compression settings, payload identity, inventory, and output checksum verification. Recorded input paths are provenance only and do not prevent resume when identical content is staged at a different path. Without `--sample-id`, the input basename determines output filenames; use a stable `--sample-id` when renamed input or GFF files must resume the same run. Resume does not continue a partially written stage.

## GFF3 mapping scope

The current mapper resolves a protein FASTA identifier against GFF3 `ID`, `protein_id`, `locus_tag`, or `Parent`. Missing and ambiguous mappings are errors. Single-CDS bacterial proteins are supported; multipart spliced CDS models are not yet supported.

FAA runs without GFF3 provide protein-relative coordinates and record that genomic coordinates are unavailable. Their short-ORF recall remains limited by the upstream annotation or gene-calling workflow.

## References

1. Larralde M. Pyrodigal: Python bindings and interface to Prodigal, an efficient method for gene prediction in prokaryotes. *Journal of Open Source Software*. 2022;7(72):4296. https://doi.org/10.21105/joss.04296
2. Hyatt D, Chen G-L, LoCascio PF, Land ML, Larimer FW, Hauser LJ. Prodigal: prokaryotic gene recognition and translation initiation site identification. *BMC Bioinformatics*. 2010;11:119. https://doi.org/10.1186/1471-2105-11-119
3. Keil B. *Specificity of Proteolysis*. Springer Berlin Heidelberg; 1992. https://doi.org/10.1007/978-3-642-48380-6
4. Gasteiger E, Hoogland C, Gattiker A, Duvaud S, Wilkins MR, Appel RD, Bairoch A. Protein Identification and Analysis Tools on the ExPASy Server. In: *The Proteomics Protocols Handbook*. Humana Press; 2005:571-607. https://doi.org/10.1385/1-59259-890-0:571
5. SIB Swiss Institute of Bioinformatics. PeptideCutter: cleavage specificities of selected enzymes and chemicals. https://web.expasy.org/peptide_cutter/peptidecutter_enzymes.html
6. Gross E, Witkop B. Selective cleavage of the methionyl peptide bonds in ribonuclease with cyanogen bromide. *Journal of the American Chemical Society*. 1961;83(6):1510-1511. https://doi.org/10.1021/ja01467a052

## License

Copyright (C) 2026 Jun Won Lee.

GenoToPep is licensed under GPL-3.0-or-later. Pyrodigal is a GPL-3.0-or-later dependency. License compatibility should be reviewed before public redistribution; this statement is not legal advice.
