How to add a covalent crosslink by hand#

Goal#

Tell tLeap to form a covalent bond between two specific atoms during the build - useful for disulfide bridges with non-standard cysteines, head-to-tail cyclic peptides, isopeptides, drug-protein covalent adducts, and any other crosslink that detectNonStandardResidues() can’t infer on its own.

Note

You usually don’t need to write custombonds by hand. parameterizeFromSpecs() walks the prepared molecule’s mol.bonds graph and emits inter-residue crosslink entries as out.custombonds for you. Reach for this recipe only when the bond isn’t in mol.bonds at all (e.g. you’re modelling a covalent adduct that wasn’t crystallised), or when you want to add bonds on top of what the canonical NCAA flow found.

Minimal example#

from htmd.builder import amber

custombonds = [
    # (atom selection for endpoint 1, atom selection for endpoint 2)
    ('segid "P0" and resid 1 and name "CA"',
     'segid "P0" and resid 50 and name "CA"'),
]

amber.build(
    mol,
    outdir="./build",
    custombonds=custombonds,
)

Each entry in custombonds is a tuple of two atom-selection strings (each must resolve to exactly one atom). At build time, amber.build translates them into tLeap bond directives.

Parameters that matter#

Parameter	What it does
`custombonds`	List of `(sel1, sel2)` atom-selection-string pairs. Each pair adds one `bond mol.<pos1>.<name1> mol.<pos2>.<name2>` to the generated `tleap.in`.
`disulfide`	Separate list of disulfide bonds with the same shape. `amber.build` auto-detects S-S within 3 Å on any pair of `CY*`-resname residues’ `SG` atoms - you only need this for non-standard cases.

Common variations#

Head-to-tail cyclic peptide#

custombonds = [
    ('segid "PEP" and resid 1 and name "N"',
     'segid "PEP" and resid 10 and name "C"'),
]
amber.build(mol, outdir="./build",
            custombonds=custombonds,
            caps={"PEP": ("none", "none")})   # suppress ACE/NME caps on the cyclised segment

Isopeptide crosslink (Gln side-chain to Lys ε-NH)#

custombonds = [
    ('segid "P0" and resid 18 and name "CD"',   # Gln18 side-chain γ-carbon
     'segid "P0" and resid 42 and name "NZ"'),  # Lys42 ε-nitrogen
]

You’ll usually also need to template both endpoints under different resnames so each side carries the right partial charges - see Build a stapled peptide for the full flow.

Multi-bond entry from `parameterizeFromSpecs`#

out = parameterizeFromSpecs(specs, prepared, outdir="./params")
amber.build(
    prepared,
    outdir="./build",
    custombonds=out.custombonds,           # list of tuples emitted by the spec inference
    topo=out.topo_paths,
    param=out.frcmod_paths,
)

When you use the canonical NCAA flow, out.custombonds already contains the right selections for every inter-residue bond detect found. You only hand-write custombonds for cases detect doesn’t cover.

Gotchas#

Each selection in a custombonds tuple must resolve to exactly one atom. If a selection matches zero or several atoms, amber.build raises a clear error. Use mol.atomselect(sel).sum() to debug.
The atoms named in a bond must already exist on the two residues - tLeap won’t add atoms, only bonds. If your endpoint is a side-chain atom that’s normally stripped (e.g. an OXT of a mid-chain residue), template the residue first so the right atoms are present.
For disulfides between residues amber.build already auto-detects (detectDisulfideBonds() matches any CY* resname with an SG atom within 3 Å), you don’t need a custombonds entry - it’ll just add a duplicate bond directive. Pass disulfide=[...] only when you want to override the auto-detection.
For a head-to-tail cyclic peptide whose terminal N-C distance is already short (< 1.35 Å in the input coordinates), amber.build’s cyclic-segment detector picks up the closure on its own and emits the appropriate bond directive without needing a custombonds entry. The minimal example above is only needed when the cyclising bond is implied by chemistry but not present in the input geometry.