Iamkaant at 22:07, 8 July 2026

2026-07-08T22:07:07Z

← Older revision		Revision as of 22:07, 8 July 2026
Line 3:		Line 3:
	'''prepare_protein.py''' is a Python 3 script that automates the receptor and		'''prepare_protein.py''' is a Python 3 script that automates the receptor and
	ligand preparation procedure that is otherwise done by hand in Chimera + vim		ligand preparation procedure that is otherwise done by hand in Chimera + vim
	before running [[Blastermaster]]. Starting from the ~~Chimera~~-minimized files		before running [[Blastermaster]]. Starting from the Maestro-minimized files
	(<code>rec_minimized_final.pdb</code> and <code>xtal_minimized_final.pdb</code>),		(<code>rec_minimized_final.pdb</code> and <code>xtal_minimized_final.pdb</code>),
	it produces the cleaned, charged, correctly-named PDB files and the		it produces the cleaned, charged, correctly-named PDB files and the
Line 29:		Line 29:
	== Quick start ==		== Quick start ==

	~~Run~~ in the directory that holds the two input files:		On Gimel, run in the directory that holds the two input files:

	<syntaxhighlight lang="bash">		<syntaxhighlight lang="bash">
	python3 prepare_protein.py		python3 /mnt/nfs/exa/work/ak87/UCSF/SCRIPTS/DOCKING/prepare_protein.py
	</syntaxhighlight>		</syntaxhighlight>

Line 38:		Line 38:

	<syntaxhighlight lang="bash">		<syntaxhighlight lang="bash">
	python3 prepare_protein.py \		python3 /mnt/nfs/exa/work/ak87/UCSF/SCRIPTS/DOCKING/prepare_protein.py \
	--rec rec_minimized_final.pdb \		--rec rec_minimized_final.pdb \
	--xtal xtal_minimized_final.pdb \		--xtal xtal_minimized_final.pdb \

Iamkaant: first commit

2026-07-08T22:04:41Z

first commit

New page

= Protein preparation for docking (prepare_protein.py) =

'''prepare_protein.py''' is a Python 3 script that automates the receptor and
ligand preparation procedure that is otherwise done by hand in Chimera + vim
before running [[Blastermaster]]. Starting from the Chimera-minimized files
(<code>rec_minimized_final.pdb</code> and <code>xtal_minimized_final.pdb</code>),
it produces the cleaned, charged, correctly-named PDB files and the
<code>working/</code> directory that Blastermaster expects.

It folds in the logic of the two legacy helper scripts
(<code>replace_his_with_hie_hid_hip.py</code> and
<code>0000_remove_hydrogens_from_pdb.py</code>), so those no longer need to be
run separately.

== Requirements ==

* Python 3 (no third-party packages; standard library only)
* Two input PDB files that have already been through the manual Chimera steps
(open structure, add hydrogens / protonate with reduce, check termini and the
protonation states of charged residues, then save):
** the receptor, saved as <code>rec_minimized_final.pdb</code>
** the crystallographic ligand, saved as <code>xtal_minimized_final.pdb</code>

The receptor file must still contain its hydrogens (including the polar
<code>HD1</code>/<code>HE2</code> on histidines) — these are required to assign
the HIS protonation states. The script removes hydrogens itself at the correct
stages.

== Quick start ==

Run in the directory that holds the two input files:

<syntaxhighlight lang="bash">
python3 prepare_protein.py
</syntaxhighlight>

Or specify paths explicitly:

<syntaxhighlight lang="bash">
python3 prepare_protein.py \
--rec rec_minimized_final.pdb \
--xtal xtal_minimized_final.pdb \
--output-dir .
</syntaxhighlight>

== What it does ==

=== Ligand pipeline ===

# Strips all header, <code>REMARK</code>, <code>CONECT</code> and <code>END</code> lines.
# Deletes all hydrogens (sphere matching does not use them).
# Writes <code>xtal-lig.pdb</code>.

=== Receptor pipeline ===

# '''Checks + fixes residue numbering.''' If the order of residues in the file does not match their residue numbers, the residue lines are reordered into ascending numeric order (the residue number is treated as ground truth). See [[#Residue numbering|Residue numbering]].
# '''Removes insertion-code letters.''' A residue such as <code>NMA 378A</code> is renumbered to the next integer (<code>379</code>), and every following residue in that chain is shifted up by one so nothing collides.
# '''Strips headers.''' Removes <code>HEADER</code> / <code>TITLE</code> / <code>REMARK</code> / <code>ANISOU</code> / <code>CONECT</code> / <code>TER</code> / <code>END</code> records, keeping only atoms.
# '''"del HC".''' Deletes hydrogens bonded to carbon (the equivalent of Chimera's <code>del HC</code>); polar hydrogens on N/O/S are kept. See [[#Removing carbon-bound hydrogens|del HC]].
# '''Caps.''' Converts <code>ACE</code> and <code>NMA</code> cap records from <code>HETATM</code> to <code>ATOM</code>, and renames the NMA cap <code>CA</code> atom to <code>CM</code>.
# '''Backbone amide H.''' Renames the backbone amide hydrogen (<code>H1</code>) of the residue following an ACE cap to <code>H</code>.
# '''Disulfides.''' Detects disulfide-bonded cysteines and renames them <code>CYS → CYX</code>. See [[#Disulfide detection|Disulfide detection]].
# Writes <code>rec_noHC.pdb</code>.
# '''HIS protonation.''' Assigns <code>HID</code> / <code>HIE</code> / <code>HIP</code> from the reduce hydrogens (<code>HD1</code> → delta, <code>HE2</code> → epsilon, both → HIP). Writes <code>rec_noHC.crg.pdb</code> (this is <code>rec.crg.pdb</code>) and copies it into <code>working/</code>.
# '''Removes all remaining hydrogens''', writing <code>rec_noH.pdb</code> (this is <code>rec.pdb</code>).

== Output layout ==

Inside <code>--output-dir</code> (default = current directory):

<pre>
xtal-lig.pdb # ligand, heavy atoms only
rec_noHC.pdb # receptor, cleaned, carbon-H removed, caps/CYX/numbering fixed
rec_noHC.crg.pdb # + HIS protonation states assigned (== rec.crg.pdb)
rec_noH.pdb # + all hydrogens removed (== rec.pdb)
rec.pdb # final no-hydrogen receptor
working/
rec.crg.pdb # the ONLY file placed in working/ (charged receptor for Blastermaster)
</pre>

The two files you normally hand to Blastermaster are
<code>working/rec.crg.pdb</code> (charged) and <code>rec.pdb</code> (no
hydrogens), together with <code>xtal-lig.pdb</code>.

== Command-line options ==

{| class="wikitable"
! Option !! Default !! Description
|-
| <code>--rec FILE</code> || <code>rec_minimized_final.pdb</code> || Chimera-minimized receptor PDB.
|-
| <code>--xtal FILE</code> || <code>xtal_minimized_final.pdb</code> || Chimera-minimized ligand PDB.
|-
| <code>--output-dir DIR</code> || <code>.</code> (current dir) || Directory to write prepared files into.
|-
| <code>--cyx "LIST"</code> || (none) || CYS residues to force to CYX, e.g. <code>"32 56 297 338"</code> or chain-qualified <code>"A32 B56"</code>. Combined with automatic detection unless <code>--no-auto-disulfide</code> is given.
|-
| <code>--no-auto-disulfide</code> || off || Disable automatic SG–SG disulfide detection; use only the residues supplied via <code>--cyx</code>.
|-
| <code>--ss-cutoff Å</code> || <code>2.5</code> || SG–SG distance (Å) below which a disulfide is called.
|-
| <code>--keep-het "LIST"</code> || (none) || <code>HETATM</code> residue names to keep, e.g. <code>"CA ZN NAG"</code>. Waters and unlisted hetero groups are dropped. ACE/NMA caps are always kept and converted to ATOM.
|-
| <code>--keep-carbon-h</code> || off || Keep hydrogens bonded to carbon (skip "del HC"). Polar hydrogens are always kept regardless.
|-
| <code>--relabel-chains</code> || off || Relabel chain IDs to sequential letters A, B, C… in order of appearance (fixes gaps such as A, C → A, B).
|-
| <code>--unify-chain-id ID</code> || (none) || Give every chain the same ID in the final <code>rec.pdb</code> only; <code>rec.crg.pdb</code> keeps proper chains for Blastermaster.
|-
| <code>--no-reorder</code> || off || Do not reorder residues by number. Insertion-code letters are still stripped.
|}

== Feature details ==

=== Residue numbering ===

Sometimes the order in which residues appear in a file does not match their
residue numbers (for example a residue numbered 14 physically appearing after
residue 378). By default the script detects this and '''reorders the residue
lines into ascending numeric order within each chain''', treating the residue
number as ground truth.

* Numbering gaps (missing residues) are preserved — they are legitimate.
* Insertion codes are then removed: <code>378A</code> becomes <code>379</code>, and every following residue is shifted up by one to avoid a collision (a cascade). Chains are handled independently.
* Use <code>--no-reorder</code> to leave the physical order untouched (letters are still stripped).

The run log reports exactly what moved and what was renumbered, e.g.:

<pre>
[check] residue numbering does NOT match sequence order in rec_minimized_final.pdb:
* chain A: LEU 14 (line 8986) is out of order -- its number is not greater than the preceding residue NMA 378A (line 8980)
-> residues will be reordered by number (see below)
reordered to numeric position: A/LEU 14
renumbered (insertion codes removed / following residues shifted):
A/NMA 378A -> 379
</pre>

=== Disulfide detection ===

Disulfides are found geometrically: the script measures the distance between
every pair of cysteine <code>SG</code> atoms (across all chains) and renames both
residues <code>CYS → CYX</code> when the distance is below the cutoff
(default 2.5 Å). It does '''not''' rely on <code>CONECT</code> records, so it works
even when connectivity is missing or incomplete.

For every cysteine it also prints the nearest-neighbour SG distance, so
borderline cases are easy to inspect:

<pre>
disulfides (auto, SG-SG < 2.5 A):
CYX A/32 -- CYX A/56 (2.05 A)
all CYS SG nearest-neighbour distances:
CYS A/32 nearest A/56 2.05 A -> CYX
CYS A/297 nearest A/338 3.44 A
</pre>

To override the geometry (e.g. to force a bridge, or to name residues yourself):

<syntaxhighlight lang="bash">
# add explicit residues to the auto-detected set
python3 prepare_protein.py --cyx "297 338"

# use ONLY the listed residues, chain-qualified, no geometry
python3 prepare_protein.py --cyx "A32 A56" --no-auto-disulfide

# loosen/tighten the distance threshold
python3 prepare_protein.py --ss-cutoff 2.3
</syntaxhighlight>

=== HIS protonation ===

Histidine protonation is read directly from the reduce-added hydrogens present
in the input:

{| class="wikitable"
! Hydrogen present !! Assigned residue name
|-
| <code>HD1</code> only (on ND1) || <code>HID</code>
|-
| <code>HE2</code> only (on NE2) || <code>HIE</code>
|-
| both <code>HD1</code> and <code>HE2</code> || <code>HIP</code>
|}

Because these markers are on nitrogen, they survive the "del HC" step, so the
assignment is unaffected by hydrogen removal.

=== Removing carbon-bound hydrogens ("del HC") ===

For each hydrogen the script finds its bonded heavy atom (the nearest
non-hydrogen atom in the same residue) and deletes the hydrogen when that atom
is a carbon. Polar hydrogens on N/O/S — including the backbone amide H, the HIS
<code>HD1</code>/<code>HE2</code>, and cap N–H atoms — are kept. This is
geometry-based rather than name-based, so it correctly handles ambiguous names
(such as <code>HD1</code> on ND1 versus a carbon) and the renamed NMA
<code>CM</code>. Pass <code>--keep-carbon-h</code> to skip this step.

=== Chains and HETATM ===

* Renumbering, disulfide detection, caps, the H1→H rename, and HIS assignment all work per chain, so multi-chain receptors are supported.
* <code>--relabel-chains</code> renames chains to sequential letters in order of appearance (e.g. A, C → A, B).
* <code>--unify-chain-id A</code> gives every chain the same ID in the final <code>rec.pdb</code> only, while <code>rec.crg.pdb</code> keeps the proper chain names — this is the trick for getting Blastermaster to renumber a 2-chain receptor.
* Waters and other hetero groups are dropped by default. Keep specific ions or cofactors with <code>--keep-het "CA ZN"</code>; kept HETATM groups are placed after the protein and are not reordered or renumbered.

== Examples ==

<syntaxhighlight lang="bash">
# Default run in the current directory
python3 prepare_protein.py

# Keep a catalytic calcium and a zinc
python3 prepare_protein.py --keep-het "CA ZN"

# Two-chain receptor: relabel A,C -> A,B and unify IDs in rec.pdb
python3 prepare_protein.py --relabel-chains --unify-chain-id A

# Force a specific disulfide set, no geometric detection
python3 prepare_protein.py --cyx "A32 A56 A297 A338" --no-auto-disulfide

# Write outputs to a separate directory
python3 prepare_protein.py --output-dir prep_4x93
</syntaxhighlight>

== Verifying the output ==

The run log is the first sanity check. Confirm that:

* the numbering check passes (or reports the fix it applied);
* the detected disulfides match what you expect;
* the HIS states match your inspection in Chimera;
* the final lines read <code>NMA CA->CM applied: True</code> and <code>NMA still CM: True</code>.

Then open <code>rec.pdb</code> and <code>working/rec.crg.pdb</code> in a viewer
and check the termini, caps, disulfides and any kept ions before running
Blastermaster.

== See also ==

* [[Blastermaster]]
* [[DOCK 3.8]]
* [[Preparing a receptor for DOCK]]

[[Category:DOCK]]
[[Category:Receptor preparation]]
[[Category:Tutorials]]

Protein Target Preparation (before Blastermaster) - Revision history

Iamkaant at 22:07, 8 July 2026

Iamkaant: first commit