Core workflow

class phenoms.SimulationSet(pdb_files=None, resid_range=None, sub_frames=None, bond_statistics_threshold=None, output_dir=None, *, trajectories=None, topology=None, topologies=None, backbone_only=True)

Bases: object

Single or replicate MD simulations with H-bond analysis.

Simple form (unchanged): pass multi-frame PDB paths via pdb_files.

Native trajectories: pass trajectories= plus topology= / topologies=, or use from_trajectories().

• backbone_only=True (default): Baker–Hubbard backbone N–O bonds (HDX-style).

• backbone_only=False: all detected donor/acceptor H-bonds.

• resid_range: filters heatmaps / manifold plots only; detection still runs on protein.

classmethod from_trajectories(trajectories, topology=None, topologies=None, *, resid_range=None, sub_frames=None, bond_statistics_threshold=None, output_dir=None, backbone_only=True)

Build a SimulationSet from native MD trajectory + topology file(s).

run(n_jobs=None, use_rust=True, *, qc=False, mdp_files=None, skip_mdp_consistency=False, qc_fail_on_nonconverged=True, qc_last_fraction=0.2, qc_mean_tolerance=0.05, qc_slope_tolerance=0.0001)

Load all trajectories, run Baker–Hubbard detection, build pivot tables. If bond_statistics_threshold was set, compute per-replicate stats and mean ± std.

Parameters:

• n_jobs (int or None) – MDTraj parallel workers per replicate (Rust path ignores this). None = default_n_jobs() (all CPUs minus two).

• use_rust (bool) – If True, use Rust extension when available. Set False if Rust returns no bonds for your topology (MDTraj fallback).

• qc (bool) – If True, run a simple replicate QC pass (RMSD convergence and optional MDP checks).

• mdp_files (list[str] or None) – Optional list of .mdp paths to compare for parameter consistency.

• skip_mdp_consistency (bool) – If True, skip MDP checks even if mdp_files are supplied.

• qc_fail_on_nonconverged (bool) – If True, raise RuntimeError when any replicate fails RMSD convergence.

export_run_artifacts(output_dir)

Write per-replicate H-bond tables, Polars-backed occupancy summaries, pivots, and a small manifest.json under output_dir/raw_data/.

Safe to call again after run() with a different path.

get_hbond_dfs()

Return list of per-replicate H-bond DataFrames (Bond Label, Frame, …).

get_occupancy_tables()

Return list of per-replicate occupancy summary DataFrames.

get_qc_report()

Return QC report dict from latest run (or None if QC was not enabled).

get_pivot_tables()

Return list of per-replicate pivot tables (bond label x frame).

get_bond_labels_sorted()

Return sorted list of all bond labels (union across replicates).

get_plot_bond_labels_sorted()

Bond labels (sorted) to use for heatmaps and manifold plots.

get_plot_pivot_tables()

Reindex each replicate pivot table to the plot bond list.

property bond_statistics

None or dict with mean_lifetimes, std_lifetimes, mean_break_frequencies, std_break_frequencies.

plot_heatmaps(save_dir=None, use_legend=False)

Plot one heatmap per replicate (whole protein or region, depending on resid_range). Bond rows = union of bonds found in this set’s replicates. For comparison of two sets with aligned bond lists, use ComparisonSet.plot_heatmaps_both() instead.

plot_top_hbonds(threshold=0.5, save_path=None)

Bar plot of top H-bonds by average lifetime (presence fraction > threshold).

plot_bond_statistics(lifetime_title='Average Lifetime of Hydrogen Bonds with Standard Deviation (Threshold = 50%)', break_title='Frequency of Breaks of Hydrogen Bonds with Standard Deviation', save_lifetime_path=None, save_break_path=None)

Bar plots for mean lifetime and break frequency with error bars (requires bond_statistics_threshold).

plot_aggregated_heatmap(save_path=None)

Average presence of each bond across all replicates.

get_fluctuating_bonds(quantile=0.9)

Bonds with highest variance (fluctuation) across frames. For single-set highlighting.

run_pca(group_labels=None, plot=True, title='PCA')

PCA on this set’s replicates (one point per replicate). Same bond space as two-set case.

run_tsne(group_labels=None, perplexity=2, plot=True, title='t-SNE', random_state=42)

t-SNE on this set’s replicates.

run_isomap(group_labels=None, n_neighbors=5, plot=True, title='Isomap')

Isomap on this set’s replicates.

run_manifold_suite(group_labels=None, perplexity=2, n_neighbors=5, random_state=42, plot=True)

PCA, t-SNE, and Isomap in one call (single-set replicates as points). group_labels: e.g. one label per replicate; default = PDB filenames.

write_structure_bfactors(pdb_path, output_path, metric='variance')

Write PDB with B-factors = per-residue metric (for PyMOL/Chimera). metric ‘variance’: highlight fluctuating residues; ‘mean’: mean occupancy. Requires biopython. Single-set only (no comparison).

class phenoms.ComparisonSet(set_a, set_b, label_a='set_a', label_b='set_b', *, output_dir=None)

Bases: object

Two simulation sets (e.g. no ligand vs ligand): compare, plot difference, write PDB, run PCA/t-SNE/Isomap.

compare(impute_threshold=0.4, flip_difference=False, donor_aggregation='sum')

export_connectivity_graph_html(output_html_path, *, graph_mode='diff', top_k_edges=30, min_abs_edge_delta=None, directed=False, title=None)

Export an interactive HTML residue connectivity graph derived from bond-level deltas.

graph_mode:

• “diff”: edges sized by abs(Δ occupancy) and colored by sign (sign consistent with flip_difference in compare()).

• “set_a” / “set_b”: edges sized by occupancy (no sign).

Parameters:

• graph_mode (str)

• top_k_edges (int)

• min_abs_edge_delta (float | None)

• directed (bool)

• title (str | None)

Return type:

None

export_connectivity_community_graph_html(output_html_path, *, graph_mode='diff', top_k_edges=40, min_abs_edge_delta=None, directed=False, resolution=1.0, title=None, community_nodes_csv_path=None, community_summary_csv_path=None)

Export a community-aware connectivity graph and optional community CSV summaries.

graph_mode:

• “diff”: edges use delta occupancy (sign preserved).

• “set_a”/”set_b”: edges use occupancy in that set.

Parameters:

• graph_mode (str)

• top_k_edges (int)

• min_abs_edge_delta (float | None)

• directed (bool)

• resolution (float)

• title (str | None)

• community_nodes_csv_path (str | None)

• community_summary_csv_path (str | None)

Return type:

None

get_comparison_df()

get_unified_bond_labels()

Union of all bond labels from both sets, sorted (for aligned heatmaps).

plot_heatmaps_both(save_dir=None, use_legend=False)

Plot heatmaps for all replicates in both sets using a unified bond list, so the same H-bonds appear in the same order on every heatmap (missing bonds = 0).

save_comparison_csv(path)

export_comparison_artifacts(output_dir)

Write raw_data/comparison.csv (bond-level Δ table) and manifest.json under output_dir. Call after compare() (or rely on get_comparison_df() which triggers compare() with defaults).

plot_difference(diff_threshold=0.2, title=None, save_path=None, diff_threshold_mode='manual', impute_small_differences=True, show_labels=True, interpolate=True, interpolation_points=900, autocorr_kappa=1.96, return_meta=False)

write_pdb_bfactors(pdb_path, output_path, value_column='Difference_clipped')

run_pca(group_labels=None, plot=True, title='PCA')

run_tsne(group_labels=None, perplexity=2, plot=True, title='t-SNE')

run_isomap(group_labels=None, n_neighbors=5, plot=True, title='Isomap')

run_manifold_suite(group_labels=None, perplexity=2, n_neighbors=5, random_state=42, plot=True)

PCA, t-SNE, Isomap on both sets’ replicates (two-group coloring by default).

phenoms.run_backbone_hbond_analysis(pdb_files, resid_range=None, sub_frames=None, n_jobs=None, plot_heatmaps=True, bond_statistics_threshold=None, output_dir=None, *, backbone_only=True)

One-shot: load replicates, run H-bond analysis, build pivot tables, optionally plot heatmaps and/or compute bond statistics.

Defaults to backbone N–O mode (backbone_only=True). Pass backbone_only=False for all-bond detection.

phenoms.simulation_set_from_dir(input_dir, prepared_dir, *, resid_range=None, sub_frames=None, bond_statistics_threshold=None, output_dir=None, frame_dt_ps=1000.0, start_ps=None, end_ps=None, apply_imaging=True, center=True, fit=True, backbone_only=True)

Build a SimulationSet from a replicate directory or set directory.

Normalizes engine-native inputs (GROMACS/OpenMM/AMBER) to multi-frame PDBs, then returns a ready-to-.run() SimulationSet. Backbone-only is the default.

phenoms.comparison_sets_from_dirs(dir_a, dir_b, prepared_dir_a, prepared_dir_b, *, resid_range=None, sub_frames=None, bond_statistics_threshold=None, output_dir_a=None, output_dir_b=None, output_dir=None, frame_dt_ps=1000.0, start_ps=None, end_ps=None, apply_imaging=True, center=True, fit=True, label_a='set_a', label_b='set_b', backbone_only=True)

Build two SimulationSet objects from two class directories (e.g., WT and MUT).

output_dir controls where the ComparisonSet’s own default bundle goes (None = auto-timestamped, False = disabled); output_dir_a/ output_dir_b control each SimulationSet’s own bundle independently.

Returns (set_a, set_b, comparison_set).

phenoms.detect_hbonds(path, *, top=None, sub_frames=None, n_jobs=None, use_rust=True, backbone_only=True)

Detect hydrogen bonds from a PDB or native trajectory path.

Parameters:

• path (str) – Trajectory / PDB path.

• top (str or None) – Topology for native MD formats (.xtc/.nc/…).

• sub_frames (int or None) – Number of frames to process; None uses all.

• n_jobs (int or None) – Parallel jobs for MDTraj fallback; defaults to default_n_jobs().

• use_rust (bool) – If True, use Rust backend when available.

• backbone_only (bool) – If True (default), keep only backbone N–O labels. If False, return all detected donor/acceptor H-bonds.

phenoms.detect_hbonds_with_occupancy(path, *, top=None, sub_frames=None, n_jobs=None, use_rust=True, backbone_only=True, output_csv_path=None)

Detect H-bonds and return both frame-level and occupancy summaries. Optionally writes occupancy summary CSV.

phenoms.default_n_jobs()

Parallel workers for MDTraj frame loop: all CPUs minus two (minimum 1).

phenoms.default_output_root()

Root directory for PHENOMS artifacts (benchmarks, test-plot runs, exports).

Resolution order: 1. Environment variable PHENOMS_OUTPUT_DIR (expanded user path). 2. Otherwise Path.cwd() / "phenom_outputs".

Return type:

Path