Bio Applied Structural Methods

Structural Determination Methods

When to Use

• Deciding whether X-ray crystallography, cryo-EM, or NMR fits a target (protein size, flexibility, membrane context, need for a ligand Kd)
• Parsing a PDB CRYST1 record to get unit cell dimensions, space group, and crystal system
• Generating symmetry-related copies of an atom/domain from a crystallographic symmetry operator
• Evaluating whether a downloaded PDB structure is good enough to use for docking, homology modeling, or mutagenesis design (resolution, R-work/R-free)
• Interpreting a PDB header's EXPDTA, REMARK 2, and REMARK 3 records

Version Compatibility

• Biopython >=1.81 (Bio.PDB.PDBParser as an alternative to manual header parsing)
• NumPy >=1.24, Python >=3.10
• wwPDB legacy PDB format (fixed-width CRYST1/REMARK records) — stable since PDB format v2; also present, differently formatted, in mmCIF (_cell.*, _refine.* tags)

Prerequisites

• pip install biopython numpy
• Basic PDB file layout (ATOM/HETATM records) — see bio-core-protein-structure
• Fractional vs. Cartesian coordinates helps but isn't required for the symmetry example below

Parsing the Unit Cell (CRYST1) and Crystal System

Goal: extract unit cell parameters and space group from a PDB file, and classify the crystal system. Approach: CRYST1 is fixed-width, not whitespace-delimited — some space groups (e.g. P 21 21 21) contain internal spaces, so split on column positions per the wwPDB spec, not str.split().

def parse_cryst1(pdb_file: str) -> dict:
    """
    Extract unit cell parameters from the CRYST1 record of a PDB file.

    CRYST1 column layout (0-indexed, fixed-width):
      cols  6-15  a (Angstrom)
      cols 15-24  b (Angstrom)
      cols 24-33  c (Angstrom)
      cols 33-40  alpha (deg)
      cols 40-47  beta  (deg)
      cols 47-54  gamma (deg)
      cols 55-66  space group
      cols 66-70  Z (copies of the asymmetric unit per cell)
    """
    with open(pdb_file) as fh:
        for line in fh:
            if line.startswith("CRYST1"):
                return {
                    "a": float(line[6:15]), "b": float(line[15:24]), "c": float(line[24:33]),
                    "alpha": float(line[33:40]), "beta": float(line[40:47]), "gamma": float(line[47:54]),
                    "spacegroup": line[55:66].strip(),
                    "Z": int(line[66:70].strip()) if len(line) > 66 and line[66:70].strip() else None,
                }
    return None


def classify_crystal_system(a: float, b: float, c: float,
                             alpha: float, beta: float, gamma: float, tol: float = 0.5) -> str:
    """Classify the 7 crystal systems from unit cell parameters (tol in Angstrom/degrees)."""
    eq = lambda x, y: abs(x - y) < tol
    is90 = lambda x: eq(x, 90.0)
    if eq(a, b) and eq(b, c) and is90(alpha) and is90(beta) and is90(gamma):
        return "Cubic"
    if eq(a, b) and is90(alpha) and is90(beta) and eq(gamma, 120.0):
        return "Hexagonal/Trigonal"
    if eq(a, b) and is90(alpha) and is90(beta) and is90(gamma):
        return "Tetragonal"
    if is90(alpha) and is90(beta) and is90(gamma):
        return "Orthorhombic"
    if is90(alpha) and is90(gamma) and not is90(beta):
        return "Monoclinic"
    return "Triclinic"


## Synthetic monoclinic example (P 1 21 1 space group)
line = "CRYST1   40.960   18.650   22.520  90.00  90.77  90.00 P 1 21 1      2\n"
import tempfile, os
with tempfile.NamedTemporaryFile(mode="w", suffix=".pdb", delete=False) as tmp:
    tmp.write(line)
    tmp_path = tmp.name
params = parse_cryst1(tmp_path)
os.unlink(tmp_path)
print(params["spacegroup"], classify_crystal_system(
    params["a"], params["b"], params["c"], params["alpha"], params["beta"], params["gamma"]))
## P 1 21 1 Monoclinic

Applying Crystallographic Symmetry Operations

Goal: generate the coordinates of a symmetry-related copy of an atom from a space-group operator. Approach: every space-group symmetry operation is x' = R @ x + t in fractional coordinates; apply the rotation matrix and translation vector from the space-group tables (e.g. International Tables for Crystallography).

import numpy as np

def apply_symmetry_operation(coords, rotation, translation) -> np.ndarray:
    """
    Apply one crystallographic symmetry operation to fractional coordinates.

    coords      : array-like, shape (3,) — fractional coordinates of an atom
    rotation    : array-like, shape (3, 3) — rotation matrix of the operator
    translation : array-like, shape (3,) — translation vector of the operator
    """
    return np.dot(rotation, coords) + translation

## 2-fold screw axis along b, space group P 1 21 1: (-x, y + 1/2, -z)
rotation_P21 = np.array([[-1, 0, 0], [0, 1, 0], [0, 0, -1]], dtype=float)
translation_P21 = np.array([0.0, 0.5, 0.0])

original = np.array([0.12, 0.34, 0.56])
sym_copy = apply_symmetry_operation(original, rotation_P21, translation_P21)
print(f"Original: {original}  ->  Symmetry mate: {sym_copy}")
## Original: [0.12 0.34 0.56]  ->  Symmetry mate: [-0.12  0.84 -0.56]

Assessing Structure Quality from the PDB Header

Goal: decide whether a downloaded PDB entry is trustworthy enough for downstream analysis (docking, mutagenesis design, homology templates). Approach: pull EXPDTA (method), resolution, and R-work/R-free from REMARK lines; a large R_free - R_work gap (>0.05-0.10) signals overfitting regardless of nominal resolution.

def parse_pdb_header(pdb_text: str) -> dict:
    """Extract experiment type, resolution, R-factors, unit cell, and space group."""
    info = {}
    for line in pdb_text.splitlines():
        if line.startswith("EXPDTA"):
            info["experiment"] = line[10:].strip()
        elif "RESOLUTION." in line:
            parts = line.split()
            for i, p in enumerate(parts):
                if p == "ANGSTROMS." and i > 0:
                    info["resolution_A"] = float(parts[i - 1])
        elif "R VALUE" in line and "WORKING" in line:
            info["R_work"] = float(line.split()[-1])
        elif "FREE R VALUE" in line and "SET" not in line:
            info["R_free"] = float(line.split()[-1])
        elif line.startswith("CRYST1"):
            info["space_group"] = line[55:66].strip()
    return info


def evaluate_structure_quality(info: dict) -> str:
    """Return a one-line quality verdict from a parsed PDB header dict."""
    res = info.get("resolution_A")
    band = ("atomic (<1.5A)" if res and res < 1.5 else
            "good (1.5-2.5A)" if res and res < 2.5 else
            "acceptable (2.5-3.5A)" if res and res < 3.5 else
            "low, side chains uncertain" if res else "unknown resolution")
    r_work, r_free = info.get("R_work"), info.get("R_free")
    overfit = ""
    if r_work is not None and r_free is not None:
        gap = r_free - r_work
        overfit = f", R_free-R_work={gap:.3f} ({'OVERFIT WARNING' if gap > 0.10 else 'OK'})"
    return f"{info.get('experiment', 'Unknown')}: {res} A ({band}){overfit}"


header = """
EXPDTA    X-RAY DIFFRACTION
REMARK   2 RESOLUTION.    2.00 ANGSTROMS.
REMARK   3   R VALUE            (WORKING SET) : 0.196
REMARK   3   FREE R VALUE                     : 0.229
CRYST1   69.850   69.850  103.180  90.00  90.00  90.00 P 41 21 2    8
"""
print(evaluate_structure_quality(parse_pdb_header(header)))
## X-RAY DIFFRACTION: 2.0 A (good (1.5-2.5A)), R_free-R_work=0.033 (OK)

Choosing a Structural Method

def recommend_method(mw_kda: float, dynamic: bool = False, membrane: bool = False,
                      need_ligand_kd: bool = False) -> str:
    """Rule-of-thumb structural method recommendation (not a substitute for expert judgment)."""
    if need_ligand_kd and mw_kda < 30:
        return "NMR (15N-HSQC titration gives per-residue CSP and Kd)"
    if mw_kda > 150 or membrane or dynamic:
        return "Cryo-EM (large/flexible/membrane targets; 3D classification separates states)"
    if mw_kda < 40:
        return "X-ray crystallography (small proteins crystallize readily; NMR also feasible)"
    return "X-ray crystallography (default for well-behaved 40-150 kDa targets)"

print(recommend_method(450, dynamic=True))     # Cryo-EM ...
print(recommend_method(15))                     # X-ray crystallography ...

Pitfalls

• CRYST1 is fixed-width, not whitespace-split: space groups like P 21 21 21 contain internal spaces; splitting on whitespace corrupts the space-group and Z fields.
• R-work always improves with more refinement parameters: only R_free (computed from reflections withheld from refinement) is a fair quality check; trust the R_free - R_work gap, not R-work alone.
• Resolution number alone is not sufficient: check R-free, and for cryo-EM also check local resolution — global resolution can hide poorly resolved flexible regions.
• Low-resolution cryo-EM maps (>4A) cannot reliably place side chains — backbone tracing may still be trustworthy.
• NMR "structures" are ensembles (multiple MODEL records in one PDB entry) representing conformational dynamics, not competing guesses of one static structure — do not just pick MODEL 1 and discard the rest.
• mmCIF vs legacy PDB: unit cell/resolution live in _cell.*/_refine.* mmCIF tags, not CRYST1/REMARK, for structures only distributed in mmCIF.

See Also

• bio-core-protein-structure — Bio.PDB SMCRA hierarchy, RMSD, DSSP
• structural-bioinformatics — Ramachandran plots, PWM/PROSITE, broader Bio.PDB workflows
• alphafold-structure-prediction — predicted-model alternative when no experimental structure exists
• bio-applied-proteomics — mass-spectrometry side of protein characterization (not covered here)