moleculekit.tools.modelling module#

moleculekit.tools.modelling.detectBackboneBreaks(mol, sel='protein', tol=2.0)#

Report every place where the protein backbone is broken.

Two residues are compared when they are adjacent in file order and share a segid and a chain, so a chain or segment boundary is never mistaken for a break. This is the check to run on a built or prepared structure: a break that survives into it is an unmodelled gap (or one the builder capped), and no job status reports it.

Parameters:
  • mol (Molecule) – The structure to check. Only the first frame is used.

  • sel (str) – Atom selection to check. The default limits it to protein residues.

  • tol (float) – Maximum C-N distance (Angstrom) still counted as a peptide bond. A real peptide bond is ~1.33 A.

Returns:

breaks – One entry per break, in file order: {"segid", "chain", "after_resid", "after_insertion", "before_resid", "before_insertion", "distance"}. after_* identifies the residue the break follows and before_* the one it precedes (the same convention as detectSequenceGaps()). distance is None when the C or N atom itself is missing, which is a break in the sense that continuity cannot be established – see moleculekit.tools.backbone.check_backbone() to rebuild such atoms.

Return type:

list of dict

Examples

>>> from moleculekit.tools.modelling import detectBackboneBreaks
>>> breaks = detectBackboneBreaks(mol)
>>> [(b["after_resid"], b["before_resid"]) for b in breaks]
[(964, 977)]
moleculekit.tools.modelling.detectSequenceGaps(mol, sequences)#

Detect missing-residue gaps per protein chain by aligning the observed sequence to the supplied full sequence.

Parameters:
  • mol (Molecule) – The gapped structure.

  • sequences (dict) – {chain: full_sequence} (e.g. from resolveFullSequences).

Returns:

  • gaps (list of dict) – Each {"chain", "after_resid", "before_resid", "missing_seq", "is_terminal"}. after_resid is the observed resid immediately before the gap (None for an N-terminal gap); before_resid is the observed resid immediately after it (None for a C-terminal gap).

  • skipped_ncaa_chains (list of str) – Protein chains skipped because they contain non-canonical residues.

  • mismatches (list of dict) – Positions where a residue is present but is not the residue the reference has there – engineered mutations, natural variants, or a reference belonging to a different construct. Each {"chain", "resid", "insertion", "reference", "observed", "ref_index"}, where reference/observed are the one-letter codes and ref_index is the 0-based position in sequences[chain], so sequences[chain][ref_index] can be patched to the observed residue to model the construct as crystallised rather than the wild type. Positions where either side is X (unknown or modified residue) are not reported: they are not evidence of a substitution.

moleculekit.tools.modelling.detectSplicedClashes(mol, new_mask, cutoff=2.0, targets='not protein')#

Report steric overlaps between the newly modelled residues and other atoms.

aceboltz models from the protein template alone, so a modelled loop can be placed into an ion / water / ligand site. This checks the new_mask atoms against targets (default the non-protein components) within cutoff.

Parameters:
  • mol (Molecule) – The spliced structure.

  • new_mask (numpy.ndarray) – Boolean mask over mol marking the newly modelled atoms.

  • cutoff (float) – Heavy-atom distance (Angstrom) below which a pair is a clash.

  • targets (str) – Atom selection for the atoms the new residues must not clash with.

Returns:

clashes{"new_chain", "new_resid", "target_resname", "target_resid", "target_segid", "min_distance"}, one per clashing new-residue / target pair, empty if clean.

Return type:

list of dict

moleculekit.tools.modelling.model_gaps(mol, sequence, segid, promod_img, minimize=False, build_sidechains=True, merge_distance=4)#

Closes residue gaps in a Molecule by sequence using ProMod3. Requires a ProMod3 Singularity image; see Notes.

This method will also mutate any residues in the Molecule that do not match the input sequence.

Parameters:
  • mol (Molecule) – The molecule containing the segment to model.

  • sequence (str) – The sequence to model.

  • segid (str) – The segment ID of the segment to model.

  • promod_img (str) – The path to the ProMod3 apptainer/singularity image. Follow the instructions at https://openstructure.org/promod3/3.4/container/singularity/ to obtain this image.

  • minimize (bool) – Whether to minimize the model after building it.

  • build_sidechains (bool) – Whether to build sidechains after building the model.

  • merge_distance (float) – The distance to merge fragments at.

Returns:

modeled_segment – The modeled segment.

Return type:

Molecule

Notes

This function requires a ProMod3 Singularity / Apptainer image. Follow the instructions at https://openstructure.org/promod3/ to obtain the image, then pass the path to the downloaded .sif file as promod_img. The function executes the modelling script inside the container via singularity exec, so Singularity or Apptainer must be available on $PATH.

Examples

>>> from moleculekit.molecule import Molecule
>>> from moleculekit.tools.modelling import model_gaps
>>> mol = Molecule("5VQ6")
>>> sequence = "HMTEYKLVVVGAGGVGKSALTIQLIQNHFVDEYDPTIEDSYRKQVVIDGETCLLDILDTAGQEEYSAMRDQYMRTGEGFLCVFAINNTKSFEDIHHYREQIKRVKDSEDVPMVLVGNKSDLPSRTVDTKQAQDLARSYGIPFIETSAKTRQGVDDAFYTLVREIRKHKEK"
>>> res = model_gaps(mol, sequence, "0", "./promod.img")
moleculekit.tools.modelling.prepareGapModellingInput(mol, sequences, gaps, outdir)#

Write the FASTA + gapped template PDB that aceboltz gapmodel consumes.

Parameters:
  • mol (Molecule) – The gapped structure.

  • sequences (dict) – {chain: full_sequence}.

  • gaps (list of dict) – The subset of detectSequenceGaps gaps the user chose to model.

  • outdir (str) – Directory to write input.fasta and template.pdb into.

Returns:

  • fasta_path (str)

  • template_path (str)

  • chain_map (dict) – {predicted_chain_label: original_chain}. aceboltz remaps template chains to "0","1",... in the order written here.

moleculekit.tools.modelling.spliceMissingResidues(mol, donor, chain_map=None, graft_flanks=1, min_identity=0.95, gaps=None)#

Insert only the newly added residues from donor into mol.

All original atoms (protein + ligands/metals/cofactors) are kept at their deposited coordinates AND with their deposited residue numbering, except the graft_flanks residues on each side of a filled gap (see below). Each modelled chain is rebuilt as its original residues plus the residues present in donor but absent from the original, including any lying beyond either terminus, so a donor observed over a wider range extends the chain; pass gaps to choose which of those runs are filled. Each inserted residue is numbered to fall between its flanking original residues, using the natural integer gap where there is room and insertion codes otherwise. A run too long for both (more residues than the free insertion codes can number) instead shifts the following residues of that chain up to open integer room, reported with a warning; residues before the run and all other chains keep their deposited numbering.

Each original protein chain is paired with the best-matching donor protein chain by sequence identity (labels are ignored), so donor chain/segid relabelling (e.g. aceboltz’s model.pdb) does not need a chain_map.

Where the missing residues go is read off the original’s backbone: a run is only inserted where the deposited chain is actually broken, never where it is covalently continuous. Two holes a few residues apart (common once incomplete residues have been dropped) would otherwise be merged into one by the sequence alignment, which pairs the residues between them with the wrong donor residues.

Parameters:
  • mol (Molecule) – The original (gapped) structure.

  • donor (Molecule) – The structure supplying the missing residues (a different crystal form, a higher-resolution entry, or an aceboltz model.pdb).

  • chain_map (dict, optional) – {donor_label: original_chain} overrides for chains that auto-pairing would get wrong or leave unpaired; the donor label is resolved against donor chain ids, then segids. Auto-pairing fills in every chain not pinned here. Default None (pure auto-pairing).

  • min_identity (float) – Minimum sequence identity (fraction of the original chain’s observed residues that match the donor) required to auto-pair an original chain to a donor chain. Original chains with no donor chain meeting this threshold are left unmodelled and reported via a warning.

  • graft_flanks (int) – Number of original residues on each side of an inserted run to also take from donor (default 1). A loop modeller rebuilds the backbone of the residues immediately flanking a gap so the new segment closes; keeping the original flanking residue instead would leave a stretched junction peptide bond, which downstream preparation reads as a chain break and caps. Grafting the flanking residues from donor (keeping their original numbering) restores a continuous backbone. Set to 0 to keep every original residue exactly.

  • gaps (list of dict, optional) – Restrict insertion to the gaps named here. Each entry is read for chain, after_resid and before_resid, the convention shared by detectBackboneBreaks() and detectSequenceGaps(), so the output of either can be filtered and handed straight in. after_resid is None for a run at the chain’s N-terminus and before_resid is None for one at its C-terminus, which is why gaps=detectBackboneBreaks(mol) fills exactly the breaks the structure has and leaves both termini as deposited. An empty list inserts nothing. Insertion codes are not matched on. Entries matching no insertable run are reported with a warning. Passing None fills every gap, terminal runs included. Call detectBackboneBreaks() before segmenting the structure: it never reports a break across a segid boundary, and segmenting starts a new segment at each break, so running it after leaves it reporting no breaks at all. detectSequenceGaps() skips a chain containing non-canonical residues entirely, so that chain contributes no entries here and is left unfilled.

Returns:

  • spliced (Molecule) – mol with the modelled residues inserted.

  • new_mask (numpy.ndarray) – Boolean mask over spliced marking the inserted atoms.

Examples

>>> from moleculekit.tools.modelling import (
...     detectBackboneBreaks, spliceMissingResidues)
>>> breaks = detectBackboneBreaks(mol)
>>> spliced, new = spliceMissingResidues(mol, donor, gaps=breaks)
>>> detectBackboneBreaks(spliced)
[]