Lesson 7 of 13 · 11 min

Measuring Distances, Angles, and Dihedrals

Interatomic distances, bond angles, and torsion angles are the raw geometric quantities behind structure validation, rotamer analysis, and Ramachandran plots. Bio.PDB exposes atomic coordinates both as NumPy arrays (Atom.coord) and as Vector objects (Atom.get_vector()), so you can compute these measurements directly or with the calc_angle and calc_dihedral helpers. This lesson works through each measurement on ubiquitin (PDB 1UBQ) and assembles a phi/psi table ready to plot.

Distances two ways

python

import numpy as np
from Bio.PDB import PDBList, PDBParser

# Fetch and parse ubiquitin; retrieve_pdb_file returns the local path pdb1ubq.ent
PDBList().retrieve_pdb_file("1UBQ", pdir=".", file_format="pdb")
parser = PDBParser(QUIET=True)
structure = parser.get_structure("1UBQ", "pdb1ubq.ent")
chain = structure[0]["A"]

ca10 = chain[10]["CA"]
ca11 = chain[11]["CA"]

# 1) Straight from coordinate vectors with NumPy
d_numpy = np.linalg.norm(ca10.coord - ca11.coord)

# 2) Bio.PDB overloads Atom.__sub__ to return the distance in angstroms
d_biopython = ca10 - ca11

print(f"{d_numpy:.3f} A (numpy)")
print(f"{d_biopython:.3f} A (Atom subtraction)")

3.768 A (numpy)
3.768 A (Atom subtraction)

Bond angles

python

from Bio.PDB.vectors import calc_angle
import numpy as np

res = chain[10]
n  = res["N"].get_vector()
ca = res["CA"].get_vector()
c  = res["C"].get_vector()

# calc_angle(v1, v2, v3) returns the angle at the middle vector v2 (the vertex).
# N-CA-C is the backbone bond angle tau; residue 10 is Gly, whose tau (~116 deg)
# runs wider than the ~111 deg backbone average, so don't read it as a typical value.
tau = calc_angle(n, ca, c)
print(f"N-CA-C angle: {tau:.4f} rad = {np.degrees(tau):.2f} deg")

N-CA-C angle: 2.0273 rad = 116.16 deg

Both calc_angle and calc_dihedral return radians, and calc_dihedral is signed on [-pi, pi]; wrap results in numpy.degrees for reporting. They consume Vector objects from Atom.get_vector(), not raw coord arrays. If you want degrees directly, call chain.atom_to_internal_coordinates() once and read residue.internal_coord.get_angle('phi' | 'psi' | 'chi1'), which returns degrees and None at chain termini.

Ramachandran torsions

python

from Bio.PDB.vectors import calc_dihedral
import numpy as np

# phi = C(i-1)-N(i)-CA(i)-C(i);  psi = N(i)-CA(i)-C(i)-N(i+1)
residues = [r for r in chain if r.id[0] == " "]  # standard residues only

rows = []
for i, r in enumerate(residues):
    phi = psi = None
    n  = r["N"].get_vector()
    ca = r["CA"].get_vector()
    c  = r["C"].get_vector()
    if i > 0:
        try:
            c_prev = residues[i - 1]["C"].get_vector()
            phi = np.degrees(calc_dihedral(c_prev, n, ca, c))
        except KeyError:
            pass
    if i < len(residues) - 1:
        try:
            n_next = residues[i + 1]["N"].get_vector()
            psi = np.degrees(calc_dihedral(n, ca, c, n_next))
        except KeyError:
            pass
    rows.append((r.id[1], r.resname, phi, psi))

for num, name, phi, psi in rows[:5]:
    phi_s = f"{phi:7.2f}" if phi is not None else "   None"
    psi_s = f"{psi:7.2f}" if psi is not None else "   None"
    print(f"{num:3d} {name}  phi={phi_s}  psi={psi_s}")

  1 MET  phi=   None  psi= 149.63
  2 GLN  phi= -91.02  psi= 138.26
  3 ILE  phi=-131.10  psi= 163.05
  4 PHE  phi=-115.99  psi= 140.23
  5 VAL  phi=-118.03  psi= 114.22

Try it yourself: Extend the loop to also compute chi1 (N-CA-CB-CG for most residues, N-CA-CB-OG1 for Thr, N-CA-CB-CG1 for Ile and Val), skipping Ala and Gly which lack a rotatable chi1. Then load the phi/psi columns into a pandas DataFrame and scatter-plot psi versus phi to render a Ramachandran diagram, and confirm that Gly residues populate regions forbidden for the other amino acids.

Programmatic Access with Biopython Bio.PDB Secondary Structure Assignment with DSSP