Skip to content

AI_CODE_GUIDE.md

Latest commit

 

History

History
448 lines (338 loc) · 23.8 KB

AI_CODE_GUIDE.md

File metadata and controls

448 lines (338 loc) · 23.8 KB
Raw

deTELpy code guide for future AI editors

This guide was written after reading the upstream deTELpy code in this checkout. It is meant to help a later AI quickly understand and safely edit the code before working on PIRATE.

Repo root: /Users/gini/.openclaw/workspace/ms-proteomics/detelpy

What deTELpy does

deTELpy has three user-facing workflows exposed through python -m deTEL <mode>:

  • rTEL: runs an open-search mass spectrometry workflow around FragPipe/MSFragger/Philosopher/Crystal-C/PTM-Shepherd/IonQuant. It produces psm.tsv and related search outputs.
  • eTEL: reads a coding-sequence FASTA and a FragPipe/Philosopher psm.tsv, detects likely amino-acid substitutions from delta masses, and outputs substitution/codon/protein count CSVs. Optionally generates plots and HTML reports.
  • mTEL: reads eTEL output CSVs plus tRNA abundance counts, fits a multinomial translation-error model with MCMC, and writes posterior traces, fitted substitution probabilities, and fitted energies.

The codebase is not heavily abstracted; most workflow logic is in a few large procedural modules. Biological constants (codon tables, amino acid lists, I/L handling) are duplicated across eTEL and mTEL, so edits there must be synchronized carefully.

Entrypoints and package shape

deTEL/__main__.py

Dispatches based on the first CLI argument:

  • no argument: imports and launches deTEL.gui.GUIRunner().main()
  • mTEL: calls deTEL.mTEL.mTEL.run(sys.argv[2:])
  • eTEL: instantiates deTEL.eTEL.workflow.eTEL.ETelRunner().run(sys.argv[2:])
  • rTEL: instantiates deTEL.rTEL.open_search.OpenSearchRunner().run(sys.argv[2:])

There is no normal --help at the top level; python -m deTEL --help is treated as an unknown mode. Use mode-specific help: python -m deTEL eTEL -h, etc.

deTEL/gui.py

Defines a Gooey GUI wrapper around the same three command-line modes. It mostly builds grouped GUI arguments and dispatches to the same runners. Treat it as a presentation layer; workflow behavior is in the mode modules.

Shared utilities

  • deTEL/utils.py: argparse helpers is_valid_file() and prep_folder().
  • deTEL/exceptions.py: currently mainly WrongSequenceTypeError for rTEL FASTA type failures.
  • deTEL/eTEL/__init__.py: CsvFileOutputColumnNames, the canonical output column-name enum used by eTEL/mTEL/reporting.
  • deTEL/eTEL/workflow/__init__.py: codon table/report constants and expected output schema constants.

Installation/dependency notes from local setup

Verified environment in this checkout:

  • Python 3.9.6 venv: .venv
  • deTELpy 0.1.13 installed editable/from source
  • Homebrew system deps: mono, openjdk@17
  • Runtime vars used:
    • PATH="/opt/homebrew/opt/openjdk@17/bin:$PATH"
    • JAVA_HOME="/opt/homebrew/opt/openjdk@17/libexec/openjdk.jdk/Contents/Home"
    • MONO_GAC_PREFIX="/opt/homebrew"

Dependency caveat: upstream spectrum_utils[iplot]==0.3.5 allowed pyteomics==5.0, which fails on Python 3.9. Local fix:

python -m pip install 'pyteomics<5'

The working environment has pyteomics 4.7.5.

See SETUP_NOTES.md and detelpy-venv-freeze.txt for exact local setup.

Data products and file contracts

rTEL input

  • Peptide FASTA database for MSFragger open search.
  • Folder containing Thermo .raw files.
  • FragPipe installation with the expected tools under tools/.

rTEL output

Within the raw-file folder and then moved into the rTEL output folder:

  • psm.tsv: Philosopher report table consumed by eTEL.
  • per-RAW .tsv and pepXML/mzML/other open-search files.
  • PTM-Shepherd outputs.
  • Optional IonQuant outputs if --ionquant is enabled.
  • Step marker files in the output folder: step-<n>-<name>-lock during execution, then step-<n>-<name>-success.

eTEL input

  • -f: coding-sequence FASTA (CDS nucleotide sequences), matching protein IDs used in the open search.
  • -psm: psm.tsv from FragPipe/Philosopher.
  • -s: rTEL/open-search output folder.
  • -o: eTEL output folder.
  • Optional -r: raw-file folder for report generation. If omitted, eTEL assumes the raw folder is the parent of -s.

The psm.tsv must have columns used by prepare() and downstream logic, including at least:

  • index spectrum IDs like rawfile.scan.scan.charge
  • Is Unique, Protein, Delta Mass, Prev AA, Next AA
  • MSFragger Localization, Protein Start, Protein End, Peptide Length
  • Peptide, Charge, Retention, Calculated Peptide Mass, PeptideProphet Probability, Intensity

eTEL output

For prefix PXD..., eTEL writes:

  • <prefix>_substitution_errors.csv
  • <prefix>_codon_counts.csv
  • <prefix>_peptide_counts.csv

substitution_errors.csv important columns:

  • protein, peptide, raw_file, modified_peptide
  • prev_aa, next_aa, charge, retention, calculated_peptide_mass
  • peptide_prophet_probability, delta_mass, protein_start, protein_end
  • localization_in_protein, codon, origin, destination, substitution, near_cognate
  • intensity, avg_base_intensity, intensity_based_error_rate

codon_counts.csv expected columns:

  • codon, base_count, error_count, detection_rate

peptide_counts.csv expected columns:

  • protein, total, erroneous, total_error_rate

mTEL input

  • Folder containing matching *_codon_counts.csv and *_substitution_errors.csv / files containing errors.
  • tRNA count CSV (Anticodon, Abundance, expected by MTEM).

mTEL output

Non-bootstrap run writes:

  • position_trace_<suffix>.csv/.pkl
  • position_<suffix>.csv/.pkl/.png
  • wobble_trace_<suffix>.csv/.pkl
  • wobble_<suffix>.csv/.pkl/.png
  • lik_trace_<suffix>.csv/.pkl/.png
  • fitted_substitution_probabilities_<suffix>.csv
  • fitted_energies_<suffix>.csv

Bootstrap run writes bootstrap mean CSV/PKL files and likelihood plots per bootstrap.

rTEL code map: open-search wrapper

Main file: deTEL/rTEL/open_search.py

Main flow: OpenSearchRunner.run()

  1. Parse/validate args.
  2. If --generate-config: copy default rTEL/configs/* to --config_dir, then stop.
  3. Resolve paths relative to the runner’s original CWD (prepare_path).
  4. Determine FASTA sequence type (determine_sequence_type) by checking sequence alphabet size. rTEL expects peptide FASTA; nucleotide FASTA causes WrongSequenceTypeError after the attempted workflow fails.
  5. Locate FragPipe root either from --fragpipe_bin_dir or by scanning PATH entries for a FragPipe bin directory.
  6. Determine tool paths under FragPipe root (determine_tool_paths), checking names and expected versions.
  7. initial_checks() verifies input folders/files and required tools/configs.
  8. Create output folder.
  9. Mutate config files in-place:
    • open_search.params: database_name, num_threads
    • crystal-c.params: fasta, thread, raw_file_location, output_location
    • shepherd.config: database, dataset
  10. Change CWD to the raw-file folder because Philosopher expects workspace-relative behavior.
  11. Run workflow steps:
    • clean Philosopher workspace
    • init workspace
    • find .raw files
    • run MSFragger (java -jar ... MSFragger ... raw_files)
    • find pepXML and .tsv outputs
    • run Crystal-C for each pepXML
    • find _c.pepXML
    • run Philosopher PeptideProphet
    • rewrite pepXML paths for uncalibrated mzML files
    • run ProteinProphet
    • run DB annotation
    • run Philosopher filtering
    • run Philosopher report -> psm.tsv
    • run PTM-Shepherd
    • clean Philosopher workspace
    • optionally run IonQuant
  12. In finally, move all non-.raw files from raw folder into output folder.

Important rTEL methods

  • run_subprocess(args): all external commands go through this; currently uses timeout=3600, captures stderr, exits with subprocess return code on failure.
  • touch_flag_files(step_name, item): decorator for workflow steps. Creates a lock file before the step and a success file after.
  • determine_tool_paths(): maps FragPipe tools/ contents to expected tool names/versions. This is a key place to edit for new FragPipe layouts.
  • step_rewrite_pepxml(): patches pepXML paths, likely fragile when file naming/layout changes.
  • step_move_all_but_raw_files(): moves outputs but skips .raw; does not overwrite existing destination files.

rTEL edit caveats

  • Upstream README says rTEL is Linux-only; on this Mac only config-generation was smoke-tested.
  • Config templates are edited in-place. If using default packaged config files directly, repeated runs may rewrite checked-in files. Prefer generating a config folder and editing that.
  • Hard-coded MSFragger JVM memory is huge (-Xmx217G) and thread values are fixed in the command apart from MSFragger config. This may need PIRATE tuning.
  • Tool/version detection assumes FragPipe directory structure and exact substrings.
  • determine_sequence_type() exits early if alphabet size <= 6; this is a heuristic and can misclassify unusual peptide FASTAs.
  • step_move_all_but_raw_files() runs even on errors, which is useful but can hide partial-output provenance unless logs/flag files are inspected.
  • step_find_raw_files() only matches lowercase .raw; uppercase .RAW files are skipped.
  • run_subprocess() has a fixed 3600-second timeout, likely too short for large searches.
  • PATH parsing is Unix-centric (: split); Windows support is only partial despite the classpath separator handling.

eTEL code map: empirical translation error detection

Main file: deTEL/eTEL/workflow/eTEL.py

Main flow: ETelRunner.run()

  1. Parse args and initialize logging.
  2. Resolve raw-file folder for optional reports.
  3. Build biological lookup tables:
    • amino-acid mass-difference dictionary (get_mass_substitution_dict)
    • codons/amino acids (utils.get_codons, get_amino_acids)
    • inverted codon table, with Leu/Ile merged for some logic.
  4. Read CDS FASTA into cds_dict with BioPython SeqIO.to_dict.
  5. Read psm.tsv as DataFrame indexed by spectrum ID.
  6. prepare() filters and annotates PSMs:
    • keep unique PSMs only
    • remove decoys and contaminants (proteins not in CDS FASTA IDs)
    • compute zero-shift peptide flag from delta-mass tolerance
    • infer protein/peptide termini
    • keep only unlocalized or singly localized modifications
    • remove zero-shift localized modifications
    • derive modified AA, modification location in peptide/protein, and raw-file basename
  7. Read danger_mods.csv and mark_danger_mods() to mark likely PTMs that overlap substitution-like masses/sites/termini.
  8. mark_substitutions() compares each PSM delta mass to all possible amino-acid substitution mass differences and marks candidate substitutions matching original modified AA and not dangerous.
  9. retain_subs_with_basepeptide() keeps only candidate substitutions where a zero-shift base peptide exists in the same protein/raw file spanning the modified codon; computes mean base intensity.
  10. Build a near-cognate mask (define_near_cognate_mask) where codons one Hamming distance away encode the destination AA, excluding synonymous codons.
  11. For each retained substitution, fetch the original codon from the CDS (fetch_substituted_codon) and split substitution into origin/destination.
  12. Select/rename/type output columns and write *_substitution_errors.csv.
  13. Calculate intensity-based error rate per substitution.
  14. Calculate protein/peptide counts and total error rate; write *_peptide_counts.csv.
  15. Count observed base codons in all filtered PSMs and substitution codons; write *_codon_counts.csv.
  16. If --generate-report, generate plots, extract spectra from RAW files, build substitution HTML report, and global report plots.

eTEL helper modules

  • utils.py: codon lists, codonify sequence helper, raw/tsv file listing.
  • extract_peptides.py: uses ms_deisotope.MSFileLoader to open Thermo .raw files, find scans by scan number, peak-pick, and return peptide dicts for spectrum rendering.
  • substitution_report.py: converts extracted peptide spectra into spectrum-utils HTML/JS renderable structures using a Jinja template.
  • dataset_report.py: plots codon counts, peptide histogram, and hyperscore distribution. Has a standalone CLI and is also called from ETelRunner.generate_report().
  • global_report.py: reads one or more eTEL output folders and CDS FASTA, computes aggregate substitution/codon summaries, and writes position histograms, heatmaps, detection-rate plots, and substitution-detected plots.

eTEL edit caveats

  • Spectrum index parsing in prepare() assumes a FragPipe-like index ending in .number.number.number; nonmatching indices will fail.
  • prepare() requires MSFragger Localization lower-case letters to mark modified residues.
  • get_mass_substitution_dict() hard-codes amino acid masses and special I/L handling. It removes L to I, I to L, and all destinations ending in L, then renames destination I as I/L. This is biologically consequential.
  • C mass is carbamidomethyl cysteine (CamCys), hard-coded.
  • The danger PTM filter uses 2 * mass_tol rather than mass_tol.
  • mark_substitutions() loops substitution dictionary keys and overwrites substitution if multiple mass differences match; order can matter if tolerance windows overlap.
  • retain_subs_with_basepeptide() requires a base peptide with Protein Start < mod_loc_in_protein < Protein End; boundaries are strict, not inclusive.
  • fetch_substituted_codon() treats mod_loc_in_protein as one-based amino-acid position and indexes codon list with mod_loc - 1.
  • near_cognate is cast to bool; if the mask contains NaN for synonymous cases, casting behavior should be checked when changing substitution semantics.
  • intensity_based_error_rate = intensity / (intensity + avg_base_intensity); if base intensity is missing/zero, downstream values may be NaN/inf.
  • Report generation needs RAW files and Mono/pythonnet/ms_deisotope; tests mock some of this.
  • Many pandas operations trigger SettingWithCopyWarning; current tests pass, but refactors should use .copy()/.loc deliberately.

mTEL code map: multinomial translation error model

Main files: deTEL/mTEL/mTEL.py, data.py, mtem.py, mcmc.py

Main flow: mTEL.run()

  1. Parse CLI args.
  2. Build output suffix; defaults to current date.
  3. List files under -f folder:
    • codon-count files are those with codon_counts in filename.
    • substitution files are those with errors in filename.
    • sorted lists are zipped and must have equal lengths.
  4. For each pair, create DataSet(substitution_file=sub, count_file=cc).
  5. If aggregate mode enabled, sum all datasets with DataSet.__add__() and fit one combined dataset.
  6. Create MTEM(anticodon_file=<tRNA count>, cell_volume=<cell vol>).
  7. Create MCMC(model=model).
  8. Initialize alpha (position penalties, length 2) and beta (wobble mismatch penalties, length 16) randomly from normal distributions.
  9. If bootstrapping, repeatedly sample datasets with replacement, run MCMC, write bootstrap mean posterior files.
  10. Else run MCMC once, then write traces/posteriors/plots, fitted substitution probabilities, and fitted energies.

DataSet

Converts eTEL outputs into the count matrix used by MCMC.

  • Reads codon counts with index_col=1, then strips unnamed columns.
  • Drops stop codons from counts.
  • Reads substitution errors and groups by codon and destination, counting rows.
  • Converts one-letter destination AA columns to three-letter names.
  • Merges I/L and L/I destinations into L.
  • Adds missing AA columns and codon rows with zeros.
  • For each codon, computes correct incorporations as base_count - substitutions_per_codon and stores them in the synonymous/correct AA column.
  • Final substitutions matrix is codon rows x amino-acid columns, integer counts.

MTEM

Biophysical/probabilistic model of tRNA competition.

  • Reads tRNA anticodon abundance CSV.
  • Builds reverse complements for codons/anticodons.
  • Drops STOP anticodons (TTA, CTA, TCA) because release-factor energies are not modeled.
  • Calculates arrival rates from tRNA abundance, diffusion coefficient, effective length, and cell volume.
  • Calculates arrival probabilities:
    • First: arrival rate over all tRNA arrival rates.
    • First_syn: arrival rate within synonymous AA group.
    • per-anticodon pairwise arrival probabilities.
  • Detection/binding model:
    • calc_dr(ac, c, position_penalty, wobble_set) = position-weighted sum of nucleotide pair penalties for codon/anticodon positions. Position 3 penalty is fixed to 1.
    • calculate_asite_binding_probability() exponentiates negative detection rates and normalizes.
    • calculate_acceptance_probability() combines binding probabilities with arrival probabilities/rates; numba-jitted inner loop.
    • calculate_all_substitution_probabilities() maps codon -> anticodon probabilities to codon -> amino-acid probabilities and optionally combines Ile/Leu.

MCMC

Metropolis sampler over MTEM parameters.

  • State:
    • alpha: log position penalties (2 params; exponentiated when used)
    • beta: wobble mismatch penalties (16 params)
    • log_likelihood
  • prepare_dataset() reorders dataset substitution matrices to match model amino-acid/codon order and removes Ile column via model THREE_LETTER_AA_NO_I.
  • Each iteration proposes alpha, evaluates log likelihood, accepts/rejects, then proposes beta, evaluates, accepts/rejects.
  • Adaptive proposal widths every 100 iterations target acceptance roughly 25-35%.
  • Trace is saved after burn-in and thinning.

mTEL edit caveats

  • aggregate = False if args.aggregate == 'n' else 'True' returns string 'True', not boolean True. It works truthily, but is sloppy.
  • CLI arg typo: num_bootsrap is misspelled throughout; keep compatibility if renaming.
  • Bootstrap mode currently calls MCMC.run() without its required logger argument and likely raises TypeError when -nb > 0.
  • wobble_<suffix>.csv is first written with posterior samples, then overwritten with posterior means in non-bootstrap mode.
  • File pairing is by sorted filename lists, not by explicit sample ID parsing. If filenames do not sort identically, datasets can be mismatched.
  • DataSet reads codon count with index_col=1, relying on current eTEL CSV layout where the first saved column is an unnamed index and second is codon.
  • DataSet.create_substitution_matrix() assumes eTEL columns codon, destination, origin.
  • DataSet and MTEM both contain codon/AA constants and I/L assumptions; keep them consistent.
  • MCMC.run() signature requires logger, but bootstrap branch in mTEL.run() calls it without logger. Bootstrap mode likely crashes unless fixed.
  • MCMC.trace = samples * [MCMCstate] initializes repeated class objects, later replaced by instances. It works after filling but is non-idiomatic.
  • plot_posterior() expects three alpha columns though current alpha has length 2; likely unused/stale.
  • MTEM.simulate_data() has simulated_codon_counts.loc[codon:, ], probably a bug (slice from current codon onward rather than current codon row).
  • calculate_incorporation_probability() has is_start, is_start = False, False and incomplete TODO logic.

Reporting code map

Dataset report (dataset_report.py)

Generates per-dataset PNGs:

  • <prefix>_codon_counts.png
  • <prefix>_peptide_counts.png
  • <prefix>_hyperscore_distribution.png

Inputs are eTEL output CSVs plus rTEL per-raw .tsv files with protein and hyperscore columns. FDR line is based on cumulative decoy ratio with hard-coded decoy prefix rev_.

Substitution report (substitution_report.py)

Takes peptide dictionaries from extract_peptides() and uses spectrum_utils to make annotated spectrum panels. It mutates peptide strings to add modification/destination context and returns renderable objects consumed by the Jinja HTML template.

Global report (global_report.py)

Aggregates across one or more eTEL output datasets in a folder:

  • Reads codonized proteome from CDS FASTA.
  • Reads all *_substitution_errors.csv files.
  • Computes substitution counts by amino acid and codon.
  • Reads codon-count CSVs into dicts by PXD ID.
  • Plots:
    • detected substitutions vs peptide counts
    • substitution position histogram
    • codon/AA heatmap
    • detection rate plot
    • substitutions detected plot

Edit risks:

  • Dataset ID parsing assumes PXD IDs in filenames/paths.
  • Several plotting functions rely on hard-coded codon/AA order and global constants.
  • read_codon_counts_as_dict() uses CSV usecols with enum-derived column names; schema changes need updates here.

GUI and entrypoint caveats

  • GUI rTEL does not expose --generate-config; use CLI for config template creation.
  • GUI rTEL uses a directory chooser for output_dir_name, although CLI treats it as a folder name under the RAW directory.
  • deTEL/__main__.py works for normal python -m deTEL <mode> usage, but its import branches are brittle if imported as a module instead of executed as __main__.
  • There is no console_scripts entry point in setup.py; assume python -m deTEL unless PIRATE adds one.

Tests and verification

Main end-to-end-ish test file:

python -m pytest deTEL/tests/test_run_deTEL.py -q

Local result: 10 passed.

Additional eTEL report tests are under deTEL/eTEL/tests/ and rTEL unit tests under deTEL/rTEL/tests/unit/.

Useful smoke commands used locally:

source .venv/bin/activate
export PATH="/opt/homebrew/opt/openjdk@17/bin:$PATH"
export JAVA_HOME="/opt/homebrew/opt/openjdk@17/libexec/openjdk.jdk/Contents/Home"
export MONO_GAC_PREFIX="/opt/homebrew"

python -m deTEL eTEL \
  -f deTEL/tests/resources/s228c_orf_cds.fasta \
  -psm deTEL/tests/resources/psm.tsv \
  -s deTEL/tests/resources/results_ionquant2 \
  -o setup_test/etel_output \
  -p PXD018591 \
  -tol 0.005

python -m deTEL mTEL \
  -f deTEL/tests/resources/results_ionquant2 \
  -r deTEL/tests/resources/tRNA_count/yeast_tRNA_count.csv \
  -o setup_test/mtel_output \
  -s 250 -p 100 -c 4.2e-17 -t 10 -b 100 -nb -1 -a n

python -m deTEL rTEL --generate-config --config_dir setup_test/rtel_config

Local results:

  • eTEL detected 11 substitutions on bundled PSM test data.
  • mTEL finished successfully with final log-likelihood around -168.33 for the smoke settings.
  • rTEL config generation created open_search.params, crystal-c.params, shepherd.config.

High-value refactor targets for PIRATE

  1. Centralize codon tables, amino-acid names, I/L handling, stop-codon handling, and mass constants into one module.
  2. Make eTEL input/output schemas explicit; validate PSM columns before running.
  3. Split ETelRunner.run() into named pipeline stages returning typed objects/dataframes.
  4. Replace rTEL in-place config mutation with copy-then-edit semantics.
  5. Make rTEL external command runner injectable/testable and log full stdout/stderr/provenance.
  6. Fix mTEL bootstrap call missing logger.
  7. Fix mTEL aggregate boolean and typo-compatible CLI naming.
  8. Pair mTEL input files by parsed dataset prefix, not sorted list position.
  9. Add deterministic random seed option for mTEL.
  10. Add JSON/YAML run manifest recording commands, input hashes, versions, and output files.
  11. Reduce hard-coded FragPipe/JVM assumptions; expose memory/tool versions as config.
  12. Add tests around non-PXD filenames if PIRATE will handle arbitrary experiment names.

Minimal mental model

  • rTEL creates psm.tsv from RAW+protein FASTA.
  • eTEL converts psm.tsv + CDS FASTA into observed substitution/codon/protein error tables.
  • mTEL converts eTEL tables + tRNA counts into fitted translation-error model parameters.
  • Reports are optional and mostly read eTEL outputs plus raw/open-search side files.

When editing, first decide which boundary changes:

  • Changing search/open-search behavior: edit rTEL.
  • Changing substitution detection/filtering: edit eTEL.
  • Changing model fitting/probabilities: edit mTEL.
  • Changing output visualizations only: edit reporting modules.
  • Changing column names or I/L semantics: update all of eTEL, mTEL DataSet, MTEM, report readers, and tests together.