🧬 GNN pLDDT Explorer¶

This tutorial shows structural biologists how to use synth-pdb's GNN quality scorer to:

1. Score structures with a single function call
2. Interpret per-residue pLDDT confidence — analogous to AlphaFold's colour scheme
3. Compare a good helix vs a distorted decoy visually
4. Run benchmark metrics (TM-score, lDDT, GDT-TS) without external tools

Prerequisites¶

pip install synth-pdb[gnn] matplotlib

---

In [ ]:

import os
import sys

IN_COLAB = 'google.colab' in sys.modules

if IN_COLAB:
    # @title Environment Setup
    !pip install -q synth-pdb biotite matplotlib numpy scipy py3Dmol openmm
else:
    sys.path.append(os.path.abspath('../../'))

print("✅ Environment configured!")

import os import sys IN_COLAB = 'google.colab' in sys.modules if IN_COLAB: # @title Environment Setup !pip install -q synth-pdb biotite matplotlib numpy scipy py3Dmol openmm else: sys.path.append(os.path.abspath('../../')) print("✅ Environment configured!")

In [ ]:

import os, sys
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.patches as mpatches
from matplotlib.colors import LinearSegmentedColormap

sys.path.insert(0, os.path.abspath('../../'))

from synth_pdb.generator import generate_pdb_content
from synth_pdb.score import score_structure, score_batch, QualityScore
from synth_pdb.benchmark_metrics import (
    extract_ca_coords, tm_score, lddt, gdt_ts, superpose_kabsch
)

print('✅ All imports successful')

import os, sys import numpy as np import matplotlib.pyplot as plt import matplotlib.patches as mpatches from matplotlib.colors import LinearSegmentedColormap sys.path.insert(0, os.path.abspath('../../')) from synth_pdb.generator import generate_pdb_content from synth_pdb.score import score_structure, score_batch, QualityScore from synth_pdb.benchmark_metrics import ( extract_ca_coords, tm_score, lddt, gdt_ts, superpose_kabsch ) print('✅ All imports successful')

Step 1 — Generate Test Structures¶

We create three structures with known properties:

Structure	Type	Expected quality
Good	Alpha-helix, 30 residues	High — ideal backbone geometry
Decoy	Random coil, 30 residues	Low — Ramachandran outliers
Mixed	Helix + beta, 30 residues	High — realistic mixed fold

In [ ]:

print('Generating structures...')
helix_pdb  = generate_pdb_content(length=30, conformation='alpha',  minimize_energy=False)
decoy_pdb  = generate_pdb_content(length=30, conformation='random', minimize_energy=False)
beta_pdb   = generate_pdb_content(length=30, conformation='beta',   minimize_energy=False)
print('✅ Done — 3 structures generated')

print('Generating structures...') helix_pdb = generate_pdb_content(length=30, conformation='alpha', minimize_energy=False) decoy_pdb = generate_pdb_content(length=30, conformation='random', minimize_energy=False) beta_pdb = generate_pdb_content(length=30, conformation='beta', minimize_energy=False) print('✅ Done — 3 structures generated')

Step 2 — Score with a Single Call¶

The score_structure() function accepts a PDB file path or an inline PDB string. It returns a QualityScore dataclass with global and per-residue information.

In [ ]:

helix_score = score_structure(helix_pdb)
decoy_score = score_structure(decoy_pdb)
beta_score  = score_structure(beta_pdb)

print('GLOBAL QUALITY SCORES')
print('='*45)
for label, r in [('Alpha-helix', helix_score), ('Random decoy', decoy_score), ('Beta-strand', beta_score)]:
    bar = '█' * int(r.global_score * 20)
    print(f'{label:14s}: {r.global_score:.3f}  {bar}  [{r.label}]')

helix_score = score_structure(helix_pdb) decoy_score = score_structure(decoy_pdb) beta_score = score_structure(beta_pdb) print('GLOBAL QUALITY SCORES') print('='*45) for label, r in [('Alpha-helix', helix_score), ('Random decoy', decoy_score), ('Beta-strand', beta_score)]: bar = '█' * int(r.global_score * 20) print(f'{label:14s}: {r.global_score:.3f} {bar} [{r.label}]')

Step 3 — Per-Residue pLDDT¶

Each residue receives a confidence score in [0, 1], colour-coded using AlphaFold's exact scheme:

Colour	Band	Score	Meaning
🔵 Blue	Very High	≥ 0.90	Backbone and side-chain likely accurate
🩵 Cyan	High	0.70–0.90	Backbone likely accurate
🟡 Yellow	Uncertain	0.50–0.70	Use with caution
🟠 Orange	Low	< 0.50	Likely disordered or incorrect

In [ ]:

PLDDT_CMAP = LinearSegmentedColormap.from_list(
    'plddt', [(0,'#FF7D45'),(0.5,'#FFDB13'),(0.7,'#65CBF3'),(0.9,'#0053D6'),(1,'#0053D6')]
)

fig, axes = plt.subplots(3, 1, figsize=(14, 7), sharex=False)
fig.suptitle('Per-Residue pLDDT Confidence (AlphaFold colour scheme)', fontsize=14, fontweight='bold')

for ax, (name, result) in zip(axes, [
    ('Alpha-helix (Good)', helix_score),
    ('Beta-strand (Good)', beta_score),
    ('Random Decoy (Bad)', decoy_score),
]):
    scores = result.per_residue
    n = len(scores)
    colors = [PLDDT_CMAP(s) for s in scores]
    ax.bar(range(1, n+1), scores, color=colors, edgecolor='none', width=1.0)
    ax.axhline(0.9, color='#0053D6', ls='--', lw=0.8, alpha=0.5, label='Very High (0.90)')
    ax.axhline(0.7, color='#65CBF3', ls='--', lw=0.8, alpha=0.5, label='High (0.70)')
    ax.axhline(0.5, color='#FFDB13', ls='--', lw=0.8, alpha=0.5, label='Uncertain (0.50)')
    ax.set_ylim(0, 1.05)
    ax.set_ylabel('pLDDT')
    ax.set_title(f'{name}  |  mean pLDDT = {np.mean(scores):.3f}  |  global = {result.global_score:.3f}')
    ax.set_xlabel('Residue number')

patches = [
    mpatches.Patch(color='#0053D6', label='Very High (≥0.90)'),
    mpatches.Patch(color='#65CBF3', label='High (0.70–0.90)'),
    mpatches.Patch(color='#FFDB13', label='Uncertain (0.50–0.70)'),
    mpatches.Patch(color='#FF7D45', label='Low (<0.50)'),
]
fig.legend(handles=patches, loc='lower center', ncol=4, fontsize=9, frameon=True)
plt.tight_layout(rect=[0,0.06,1,1])
plt.show()
print('💡 Blue = high confidence (like AlphaFold DB).  Orange = disordered or bad geometry.')

PLDDT_CMAP = LinearSegmentedColormap.from_list( 'plddt', [(0,'#FF7D45'),(0.5,'#FFDB13'),(0.7,'#65CBF3'),(0.9,'#0053D6'),(1,'#0053D6')] ) fig, axes = plt.subplots(3, 1, figsize=(14, 7), sharex=False) fig.suptitle('Per-Residue pLDDT Confidence (AlphaFold colour scheme)', fontsize=14, fontweight='bold') for ax, (name, result) in zip(axes, [ ('Alpha-helix (Good)', helix_score), ('Beta-strand (Good)', beta_score), ('Random Decoy (Bad)', decoy_score), ]): scores = result.per_residue n = len(scores) colors = [PLDDT_CMAP(s) for s in scores] ax.bar(range(1, n+1), scores, color=colors, edgecolor='none', width=1.0) ax.axhline(0.9, color='#0053D6', ls='--', lw=0.8, alpha=0.5, label='Very High (0.90)') ax.axhline(0.7, color='#65CBF3', ls='--', lw=0.8, alpha=0.5, label='High (0.70)') ax.axhline(0.5, color='#FFDB13', ls='--', lw=0.8, alpha=0.5, label='Uncertain (0.50)') ax.set_ylim(0, 1.05) ax.set_ylabel('pLDDT') ax.set_title(f'{name} | mean pLDDT = {np.mean(scores):.3f} | global = {result.global_score:.3f}') ax.set_xlabel('Residue number') patches = [ mpatches.Patch(color='#0053D6', label='Very High (≥0.90)'), mpatches.Patch(color='#65CBF3', label='High (0.70–0.90)'), mpatches.Patch(color='#FFDB13', label='Uncertain (0.50–0.70)'), mpatches.Patch(color='#FF7D45', label='Low (<0.50)'), ] fig.legend(handles=patches, loc='lower center', ncol=4, fontsize=9, frameon=True) plt.tight_layout(rect=[0,0.06,1,1]) plt.show() print('💡 Blue = high confidence (like AlphaFold DB). Orange = disordered or bad geometry.')

Step 4 — Batch Scoring¶

When scoring many structures, use score_batch() — the model is loaded once for the whole list.

In [ ]:

print('Generating 10 structures for batch scoring...')
pdbs = [
    generate_pdb_content(length=20+i*3, conformation='alpha', minimize_energy=False)
    for i in range(5)
]
pdbs += [
    generate_pdb_content(length=20+i*3, conformation='random', minimize_energy=False)
    for i in range(5)
]

results = score_batch(pdbs)
print(f'Scored {len(results)} structures\n')
print(f'  High Quality: {sum(1 for r in results if r.label=="High Quality")}')
print(f'  Low Quality:  {sum(1 for r in results if r.label=="Low Quality")}')

print('Generating 10 structures for batch scoring...') pdbs = [ generate_pdb_content(length=20+i*3, conformation='alpha', minimize_energy=False) for i in range(5) ] pdbs += [ generate_pdb_content(length=20+i*3, conformation='random', minimize_energy=False) for i in range(5) ] results = score_batch(pdbs) print(f'Scored {len(results)} structures\n') print(f' High Quality: {sum(1 for r in results if r.label=="High Quality")}') print(f' Low Quality: {sum(1 for r in results if r.label=="Low Quality")}')

Step 5 — Structural Benchmark Metrics¶

These metrics let you compare any two structures — for example, an AI prediction vs a reference. They are the same metrics used in CASP competitions and the AlphaFold2 paper.

Metric	Scale	Same fold?	Notes
TM-score	0–1	> 0.5	Length-independent, preferred for fold comparison
GDT-TS	0–1	> 0.5	CASP standard metric
lDDT	0–1	—	Per-residue; no superposition needed
Cα-RMSD	Å (lower=better)	< 2 Å	Sensitive to outliers

In [ ]:

# Use helix as 'reference' and a noisy version as 'prediction'
rng = np.random.default_rng(42)
ca_ref = extract_ca_coords(helix_pdb)
noise  = rng.normal(0, 1.5, ca_ref.shape).astype(np.float32)
ca_noisy = ca_ref + noise   # simulated imperfect prediction

# Compute all metrics
tm     = tm_score(ca_noisy, ca_ref)
gdt    = gdt_ts(ca_noisy, ca_ref)
lddt_r = lddt(ca_noisy, ca_ref)
_, rmsd = superpose_kabsch(ca_noisy, ca_ref)

print('STRUCTURAL COMPARISON METRICS')
print('='*45)
print(f'  TM-score : {tm:.4f}  ({"same fold" if tm>0.5 else "different fold"})')
print(f'  GDT-TS   : {gdt:.4f}')
print(f'  lDDT     : {np.mean(lddt_r):.4f}  (mean over {len(lddt_r)} residues)')
print(f'  Cα-RMSD  : {rmsd:.2f} Å')

# Use helix as 'reference' and a noisy version as 'prediction' rng = np.random.default_rng(42) ca_ref = extract_ca_coords(helix_pdb) noise = rng.normal(0, 1.5, ca_ref.shape).astype(np.float32) ca_noisy = ca_ref + noise # simulated imperfect prediction # Compute all metrics tm = tm_score(ca_noisy, ca_ref) gdt = gdt_ts(ca_noisy, ca_ref) lddt_r = lddt(ca_noisy, ca_ref) _, rmsd = superpose_kabsch(ca_noisy, ca_ref) print('STRUCTURAL COMPARISON METRICS') print('='*45) print(f' TM-score : {tm:.4f} ({"same fold" if tm>0.5 else "different fold"})') print(f' GDT-TS : {gdt:.4f}') print(f' lDDT : {np.mean(lddt_r):.4f} (mean over {len(lddt_r)} residues)') print(f' Cα-RMSD : {rmsd:.2f} Å')

In [ ]:

# Plot per-residue lDDT
fig, ax = plt.subplots(figsize=(12, 4))
ax.bar(range(1, len(lddt_r)+1), lddt_r, color='steelblue', edgecolor='none')
ax.axhline(np.mean(lddt_r), color='red', ls='--', lw=1.5, label=f'Mean lDDT = {np.mean(lddt_r):.3f}')
ax.set_ylim(0, 1.05)
ax.set_xlabel('Residue number')
ax.set_ylabel('lDDT score')
ax.set_title('Per-Residue lDDT: Reference vs Noisy Prediction')
ax.legend()
plt.tight_layout()
plt.show()
print('💡 Low lDDT residues have local distance geometry that differs from the reference.')

# Plot per-residue lDDT fig, ax = plt.subplots(figsize=(12, 4)) ax.bar(range(1, len(lddt_r)+1), lddt_r, color='steelblue', edgecolor='none') ax.axhline(np.mean(lddt_r), color='red', ls='--', lw=1.5, label=f'Mean lDDT = {np.mean(lddt_r):.3f}') ax.set_ylim(0, 1.05) ax.set_xlabel('Residue number') ax.set_ylabel('lDDT score') ax.set_title('Per-Residue lDDT: Reference vs Noisy Prediction') ax.legend() plt.tight_layout() plt.show() print('💡 Low lDDT residues have local distance geometry that differs from the reference.')

Step 6 — Understanding the Difference: GNN vs Ramachandran¶

The GNN scorer and Ramachandran analysis are complementary:

Metric	Captures	Speed	Per-residue?
Ramachandran %	Backbone dihedral geometry	Fast	✅ (as plot)
GNN global score	Full contact topology	< 1 ms	✅ (pLDDT)
TM-score vs ref	Fold accuracy vs a known structure	Fast	❌
lDDT vs ref	Local accuracy vs a known structure	Fast	✅

When to use each:

• Use Ramachandran % for PDB deposition (required by journals)
• Use GNN pLDDT for rapid screening of large structure sets
• Use TM-score / lDDT when you have a reference structure to compare against

In [ ]:

# Side-by-side: GNN pLDDT vs computed lDDT for the helix
# lDDT of helix vs itself = 1.0 everywhere; vs noisy = varies
plddt_vals = np.array(helix_score.per_residue)
n_align = min(len(plddt_vals), len(lddt_r))

fig, axes = plt.subplots(1, 2, figsize=(14, 4))

axes[0].bar(range(1, n_align+1), plddt_vals[:n_align], color='#0053D6', alpha=0.8)
axes[0].set_ylim(0, 1.05)
axes[0].set_title('GNN pLDDT (internal model confidence)')
axes[0].set_xlabel('Residue')
axes[0].set_ylabel('Confidence')

axes[1].bar(range(1, n_align+1), lddt_r[:n_align], color='steelblue', alpha=0.8)
axes[1].set_ylim(0, 1.05)
axes[1].set_title('Computed lDDT vs noisy reference')
axes[1].set_xlabel('Residue')
axes[1].set_ylabel('lDDT')

plt.suptitle('GNN pLDDT vs Computed lDDT', fontsize=13, fontweight='bold')
plt.tight_layout()
plt.show()
print('Note: GNN pLDDT is a model confidence score; lDDT requires a reference structure.')

# Side-by-side: GNN pLDDT vs computed lDDT for the helix # lDDT of helix vs itself = 1.0 everywhere; vs noisy = varies plddt_vals = np.array(helix_score.per_residue) n_align = min(len(plddt_vals), len(lddt_r)) fig, axes = plt.subplots(1, 2, figsize=(14, 4)) axes[0].bar(range(1, n_align+1), plddt_vals[:n_align], color='#0053D6', alpha=0.8) axes[0].set_ylim(0, 1.05) axes[0].set_title('GNN pLDDT (internal model confidence)') axes[0].set_xlabel('Residue') axes[0].set_ylabel('Confidence') axes[1].bar(range(1, n_align+1), lddt_r[:n_align], color='steelblue', alpha=0.8) axes[1].set_ylim(0, 1.05) axes[1].set_title('Computed lDDT vs noisy reference') axes[1].set_xlabel('Residue') axes[1].set_ylabel('lDDT') plt.suptitle('GNN pLDDT vs Computed lDDT', fontsize=13, fontweight='bold') plt.tight_layout() plt.show() print('Note: GNN pLDDT is a model confidence score; lDDT requires a reference structure.')

🎓 Key Takeaways¶

1. score_structure(path_or_string) — one call gives you global quality + per-residue pLDDT
2. pLDDT > 0.9 (blue) = Very High confidence — same as AlphaFold's blue regions
3. TM-score > 0.5 means two structures share the same global fold
4. lDDT is the per-residue local accuracy metric — no superposition needed
5. GNN pLDDT and computed lDDT are complementary:
   • GNN pLDDT: fast screening without a reference
   • Computed lDDT: accurate local comparison when a reference is available

Next Steps¶

• API Reference: score_structure()
• API Reference: benchmark metrics
• Scientific Background: GNN Structure Scoring
• Try the AlphaFold vs NMR Dynamics tutorial

In [ ]:

print('Tutorial complete! 🎉')
print()
print('Quick reference:')
print('  from synth_pdb.score import score_structure')
print('  result = score_structure("my_protein.pdb")')
print('  print(result.global_score, result.label)')
print('  print(result.per_residue)   # pLDDT per residue')

print('Tutorial complete! 🎉') print() print('Quick reference:') print(' from synth_pdb.score import score_structure') print(' result = score_structure("my_protein.pdb")') print(' print(result.global_score, result.label)') print(' print(result.per_residue) # pLDDT per residue')