Does this peptide actually bind PEPR1, PEPR2 and BAK1?
This is an atomic-resolution structural model of molecular recognition: given one peptide sequence and three receptor-side molecules, it predicts the complex they form, grades that prediction against a real crystal structure, and hands the result back to peptide design. It answers whether a sequence can be bound — not whether binding triggers a response. Five zones of work, each built around one or two load-bearing modules, with supporting tools feeding them and a redesign loop running back from the verdict.
Query peptide
BoPEP4 · 23 aa · net +4
Receptor
PEPR1 — PDB 5GR8
Paralogue
PEPR2 — model only
Co-receptor
BAK1 — PDB 4MN8
Evidence flag
◇ inferred, not solved
Which protectant ran which part
The dots on each card show which protectants passed through that stage, so overlaps stay visible where more than one did.
Two branches running at once, meeting only at docking: each protectant on the left, whatever it has to act on or be modelled against on the right. Only two of the nine molecules here are experimental structures of the thing they represent.
Five checks around one verdict, failing for unrelated reasons. Which one fails is the diagnosis for the next design.
What this model cannot tell you
The limits are part of the model. A well-scored, stereochemically clean complex says the sequence is capable of the molecular recognition step — the one thing this kind of modelling is able to answer. It does not say the peptide is taken up, survives the apoplast, or that binding fires the downstream immune output. Pattern-triggered immunity is read out in the cell, not in the coordinates, so the redesign loop can only ever optimise for predicted binding, never for function.
Established here
Shape complementarity, plausible interface chemistry, a ΔG estimate, and reproduction of known anchor contacts.
Binding assay (ITC / MST / SPR), receptor co-immunoprecipitation, and a functional defence-marker readout.
Two parts of the map are flagged ◇ inferred rather than solved: the PEPR2 structure, and every PEPR–BAK1 interface. Both rest on homology plus convergent conservation evidence. Report them as models, and say so in the figure legend.
Dashed slots are figures this methodology still needs; each one names what it would show.
1Introduction
1.1 Background
Computational tools give real insight into protein engineering and act as decision-making engines for wet lab work, but the process is easy to oversimplify: convenient analysis of inputs and outputs in place of thorough analysis, and no methodology tailored to the project it serves.
Their outputs depend heavily on the quality of inputs, parameters and underlying assumptions, and that dependence is worse for short, disordered peptides such as Brassica oleracea Plant Elicitor Peptide 4 (BoPEP4), one of our three protectants alongside Late Embryogenesis Abundant protein 14 (LEA14) and 1-aminocyclopropane-1-carboxylate deaminase (ACC deaminase). A workflow that holds up has to account for how and why these tools generate their predictions, not merely run them.
Our docking methodology predicts and evaluates the ternary signalling complex of the receptor-like kinase PEPR1, the co-receptor kinase BAK1, and the Brassica-derived peptide ligand BoPEP4. We designed the workflow on the principles of integrative structural biology established by Alexandre Bonvin at Utrecht University, whose lab produced most of the software below.
How to read this page
Every stage is written twice. The opening text is the method itself: what the step is for, what it takes in, and what it can and cannot tell you. The panels under it are what we got when we ran it on each protectant.
We started from the stress we needed to fix, not from a molecule we liked, and we searched gene families instead of single proteins. A family is a set of sequences that arose by duplication, so its members usually share a mechanism even when their triggers differ. Four criteria cut the shortlist: does the mechanism match the stress we are targeting, can we assay it in the wet lab we have, can we engineer it, and is it novel enough to be worth the risk. A candidate that fails one is written down as considered and excluded. The exclusions set the edges of the design space, so they belong in the record as much as the three we kept.
[Insert PROPEP family tree, high quality image]
The PEP family, and why BoPEP4 over AtPEP1
The BoPROPEP family tree identified candidate PEP peptides including BoPEP4, AtPEP1 and AtPEP3, and let us place BoPEP4 in relation to Arabidopsis rather than taking the peptide directly from the plant. PEP-family peptides are perceived by the same pair of receptor kinases, PEPR1 and PEPR2, so receptor recognition is conserved across the family: a reasonable first signal for which candidates behave similarly, before any experiment.
We did not assume the best-characterised peptide would be the best target. AtPEP1 is known mainly for biotic stress, while the project needed measurable abiotic stress, which weather data and sensors detect well. Drought is prominent in Taiwan and is directly signalled by BoPEP4. Abiotic focus also gives direct plant phenotypes such as root length and biomass, where biotic assays would add a second organism and more confounders.
Finally, BoPEP4's variable N-terminal region alongside its conserved C-terminal motif (SSG-x2-G-x2-N) allows targeted N-terminal edits while preserving the conserved core, and its relative obscurity means the sequence variants are genuinely uncharacterised.
The LEA family, and why LEA14
LEA proteins are associated with abiotic stress and can show substantial intrinsic disorder: some lack a single stable fold and adopt different conformations depending on environment or partner. We therefore separated ordered or globular members from intrinsically disordered ones, because mutation methods designed for folded proteins cannot be read the same way for disordered peptides, and the family contains both.
We chose the ordered, globular LEA14, with demonstrated protection against drought, salinity and osmotic stress. Its defined structure gives a stable scaffold for the protein-protein interactions involved in molecular stabilisation. Unlike BoPEP4, it does not depend on recognising one dedicated receptor; its stabilisation function can involve many cellular proteins, so a larger structured protein suits the question better. A disordered LEA would instead act through broader, less specific molecular-shield mechanisms.
It also gives a measurable phenotype for wet lab testing, and at over 150 residues it poses a fundamentally different design problem from a short receptor-binding peptide.
ACC deaminase
Suggested by Stanford bioengineering assistant professor Jenn Brophy in interview. It fills a real gap: both LEA14 and BoPEP4 fail against heat stress, one of Taiwan's largest abiotic stresses.
The enzyme is largely specific to ACC, the direct precursor of the plant hormone ethylene. Acting on that one signalling molecule lets it reduce stress-induced ethylene accumulation, and so the growth inhibition that follows, without interfering with other cellular chemistry. It also functions both inside plant cells and outside the roots, which widens the application routes.
Considered and set aside: trehalose, proline, CAP peptides
Trehalose was set aside despite its established role in osmotic tolerance: it is already heavily engineered in synthetic biology and offers little differentiation. Proline was excluded because excessive accumulation is associated with cytotoxic rather than purely protective effects, a risk for a system producing protectant continuously. CAP-derived peptides were deprioritised on timeline, since synthesis and validation would not fit the development window.
[potentially add graph: the four screening criteria scored for every candidate considered, including the ones excluded, as a matrix with the three survivors highlighted, purpose: to show that the exclusions define the design space as much as the selections do, which the prose claims but nothing on the page currently demonstrates]
2.2 Literature: understanding the mechanism
The reading supplies the parameters every later step needs, and it tells us whether a design is new. We used three kinds of source: sequence databases for the canonical sequence and the domain boundaries, structural databases for coordinates and for the residue range that was actually resolved, and primary functional papers for the assays and the effect sizes. With only the first two you can build a design, but you cannot say what it means.
Tang et al. 2015
Resolved the AtPEP1-receptor structure and showed the peptide binding in an extended conformation. This is why our docking inputs treat BoPEP4 as a set of possible conformational states rather than one rigid structure, and it is the reference the receptor-binding interface is read against.
Yu et al.
Covers broader strategies for designing and engineering plant-derived peptides. It supported truncation and residue substitution as exploratory axes, and we used it to identify which modifications would produce interpretable changes, since not every mutation is equally informative.
Pearce et al. 2008
An alanine scan of AtPep1(9–23) and the basis for several of our experimental design concepts. It motivated systematic peptide modification to establish which residues contribute to activity, and it later settled a design question directly: [A18] is fully active in the alkalinization assay at 2.5 nM, ten times the EC50.
Truncation series, as a method
A truncation series progressively removes portions of a sequence, so comparing variants of different lengths shows whether a region is necessary for the expected interaction. Truncation was therefore both an engineering strategy and a mechanistic probe.
Alanine scanning, as a method
Replaces selected residues with alanine, keeping a simple side chain while removing the original chemistry. The purpose is diagnostic: a substitution that degrades the measured interaction marks a load-bearing position, and one that does not marks a degenerate position available for other modification. With truncation it is the conceptual foundation of the mutation-design workflow.
[potentially add model: the PEPR1 signalling cascade from peptide binding through BAK1 recruitment to the Ca2+ / ROS / MAPK outputs, with the single step this pipeline models boxed off from the rest, purpose: to draw the boundary of the modelling claim in one picture, so a judge sees immediately that we model recognition and not response]
2.3 Interviews
An interview is only useful to us later if the record keeps three things apart: what the expert stated, what they recommended, and what we decided afterwards. They collapse into one line easily, and then nobody can tell whose judgement it was. The questions have to be narrow too. A question with no edges comes back with an answer we cannot act on.
Dr. Jenn Brophy, Stanford Bioengineering
Her interview shaped decisions beyond the choice of ACC deaminase. Asked how to apply stress-sensor and protectant-synthesis logic so plants tolerate heat and salinity, she recommended an OR-gate architecture, and advised working backward from application needs: select a logic gate that is activated under the conditions you care about rather than starting from circuit design. She pointed to two- and three-input logic gates in E. coli as a resource.
Her ACC deaminase recommendation was scoped to drought and heat. She made no direct claim about salinity, so we researched the OR-gate further to see whether it extends to both heat and salinity.
2.4 Protectant families: sequence alignment
Comparing sequences is the cheapest first pass at conservation, variability and possible engineering sites, and it costs nothing to run before any structure exists. What it gives is candidate positions, not a ranking of them. Comparison on its own cannot tell a position that varies because it does not matter from one where two peptides solve the same problem in different ways.
PEP sequence alignment and BoPEP4
Multiple AtPEP sequences were aligned to find shared residues and conserved patterns, then compared with BoPEP4. That gives both halves of the engineering brief at once: conserved residues are features to preserve, and the positions where they differ are candidates for mutation and truncation. The AtPEP1–PEPR1 structure provides the context for reading those conserved residues, especially in the C-terminal region that makes concentrated receptor contacts.
Results · BoPEP4 alignment findings
23 positions, 10 identical to AtPEP1
Position
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
BoPEP4
G
I
L
I
G
S
K
K
R
P
R
E
P
H
S
S
G
K
P
G
G
H
N
Identity
•
•
•
•
•
•
•
•
•
•
AtPEP1
A
T
K
V
K
A
K
Q
R
G
K
E
K
V
S
S
G
R
P
G
Q
H
N
Shaded columns are identical between the two peptides. All ten fall in the C-terminal half; residues 1–6 share nothing.
Both peptides are 23 residues and align without gaps, so every position has a direct counterpart and numbering transfers one to one. This is what allows AtPEP1 structural data to be read onto BoPEP4 at all.
Ten positions are identical: 7, 9, 12, 15, 16, 17, 19, 20, 22 and 23, all in the C-terminal region. Identity is 10 of 23 across the full peptide, 10 of 17 within residues 7–23, and zero across residues 1–6. That boundary matters because an independent method draws it in the same place: 5GR8 resolves AtPep1 only from residue 7, and our contact analysis found both PEPR1 and PEPR2 engage residues 7–23 while the first six make no receptor contact.
BoPEP4 also conforms to the conserved motif SSG-x2-G-x2-N at positions 15–23: Ser15, Ser16, Gly17, Gly20 and Asn23 are all present. The two positions differing inside that window, 18 and 21, both fall in the motif's degenerate positions. The peptide satisfies the family signature and differs only where the family permits variation.
The N-terminal divergence is chemical as well as positional. AtPEP1 opens with ATKVKA, charged; BoPEP4 opens with GILIGS, uncharged and hydrophobic. The whole segment differs in character, not at scattered positions.
Across the three engineering axes the alignment fixed a starting point for each. N-terminal truncation is supported by two independent lines, zero conservation and zero receptor contact across 1–6. Point substitution has a defined candidate set inside 7–23, with positions 18 and 8 differing from AtPEP1. The conserved repeat analysis has a defined unit at 15–23. Which of those candidates survived is reported in Sprints and Results.
Why we did not run a family alignment for LEA14
Sequence conservation alone was less informative for the questions being asked. The LEA family contains substantial structural diversity, including both intrinsically disordered and folded members, so similar sequence patterns do not necessarily correspond to equivalent structural or functional roles.
LEA14 was instead analysed as a specific folded protein: its experimentally characterised structure, residue environment, interaction surfaces, and the vulnerabilities created by our delivery route. That let us connect each proposed mutation to a defined structural hypothesis. The two protectants therefore needed different forms of sequence analysis.
[potentially add graph: per-position conservation across the PROPEP family from the MUSCLE alignment, plotted along the 23 residues with the ten BoPEP4/AtPep1 identities marked, purpose: to show the conserved signal is concentrated in the C-terminal half, which is the structural justification for threading and for where restraints were placed]
2.5 Engineering axes
Our substitutions split in two, and which kind one is decides how we justify it. An interface substitution changes a contact, and we take the new residue from a homologue whose structure is solved, so it is a residue already seen to work in that position. A liability substitution removes a predicted failure mode: one chemical feature that one process condition is expected to attack.
BoPEP4 axes
Truncations
The uncharged six-residue GILIGS segment carries zero conservation and zero receptor contact, making the N-terminus the only end that can be trimmed. The C-terminus is off limits: Asn23's free carboxylate is the anchor, and deleting the last residue is documented to compromise binding. Fewer residues to synthesise and export also implies faster secretion, so shortening is a manufacturability gain if binding holds. The wet lab took this up as a testable claim and is measuring secretion for truncated constructs against wild type.
Point substitutions
The axis came from a simple idea: BoPEP4 residues not in contact with PEPR1 could be swapped for the residue AtPEP1 carries at the same position, restoring residues present in the literature-backed 5GR8 structure. Alignment identified positions 18 and 8 as the candidate set inside the contact region.
Conserved repeats
BLAST scans suggested the conserved 15–23 region drives receptor binding, motivating a tandem repeat analysis. The countervailing risk is that lengthening the peptide makes it less active-site specific, filling other parts of the receptor.
What each axis returned, including the designs we withdrew and why, is in Sprints and Results.
LEA14 axes
Point substitution
These address liabilities identified in the sequence, not contacts at an interface. LEA14 contains a single cysteine, Cys58, at the start of β3 with no partner. In the reducing plant cytosol an unpaired thiol is inert, but our route sends the protein through the oxidising extracellular compartment where it can form intermolecular disulfides and dimerise. C58S substitutes a serine of almost identical size that cannot. The test is a gel run with and without reducing agent, so the hypothesis fails visibly if a paired band persists.
The sequence also contains exactly one Asp-Pro bond, at 78–79, the most acid-labile peptide bond type. D78E removes it while preserving the negative charge. This is insurance against a process condition, not a gain in activity, and was ranked as such.
A third group targets the flexible N-terminal extension at 2–19. Helical wheel analysis gives that segment a hydrophobic moment near 0.37, periodic but below the range typical of membrane-binding amphipathic helices. Face-consistent substitutions could raise it into that range and turn a flexible tail into a membrane anchor. We record this as a hypothesis, and built the N-terminal truncation as its matched negative control.
We deliberately did not pursue surface hydrophilicity substitutions. LEA14's GRAVY of +0.056 is unusual for a stress-protective protein and could be lowered, but doing so works against the fold, and the bulk hydration shielding it would serve is the mechanism our feasibility analysis argues is unavailable at any deliverable concentration.
K-segment translocation
LEA14 belongs to the atypical folded LEA_2 group and carries none of the repeat motifs characterising the disordered groups. The K-segment is the best-characterised membrane-protection module in the family and occurs in dehydrins, not here. The design fuses one or two K-segments, or Group 3 LEA 11-mer repeats, onto the folded core, giving a stable scaffold with a disordered arm.
The precedent is documented: the mitochondrial LEA protein LEAM is disordered in solution and folds into an amphipathic helix that binds and stabilises membranes as the cell dries. If the motif carries that behaviour rather than the whole protein, grafting it should carry the behaviour with it. Two constraints bound the graft: added mass reduces the fraction crossing the cell wall, and every fusion junction is a new site to audit for degron compliance and protease susceptibility.
Truncations
The published solution structure defines residues 24–145. The first 23 have no fixed conformation in the ensemble and are where a protease attacks first. Removing them does three things at once: removes the most exposed proteolytic target, lowers mass and hydrodynamic radius, and leaves the folded β-sandwich core intact. It is also the negative control for the amphipathicity hypothesis above, so the two constructs read against each other. The test is a trypsin challenge comparing how quickly full-length and truncated bands disappear.
[potentially add model: the threaded BoPEP4 pose with the engineering-axis positions coloured by category, interface substitutions in one colour and liability substitutions in another, purpose: to show the two kinds of substitution occupy different regions of the peptide and are therefore justified on different evidence]
We based the input on two sources: 5GR8, the crystal structure of PEPR1 with AtPEP1, and AlphaFold. 5GR8 confirmed the extended conformation of AtPEP1 along the inner surface of the superhelical PEPR1 LRR, supporting the inference that BoPEP4, from the same PEP family, adopts a comparable extended structure when bound. Structures were generated with AlphaFold 3, cross-checked against 5GR8, and used as docking inputs, with 5GR8 itself as the control.
Run on the free peptide, AlphaFold gives us nothing we can dock: fraction_disordered = 1.0, pTM 0.18–0.35, and five seeds returning five different secondary structures. For a peptide that only orders when it binds, that is the expected answer and not a failed run. The route we took instead is threading onto the AtPEP1 backbone of 5GR8 chain J, which the gapless alignment allows: the backbone transfers unchanged and we rebuild only the side chains at the positions that differ.
The receptor trace is PEPR1 chain A of 5GR8; the accent trace is the peptide in its groove. Threading grafts the BoPEP4 sequence onto this backbone, so the two peptides have identical Cα coordinates by construction (deviation 0.00 Å across all 17 resolved positions, 7–23). Ten positions carry the same residue in both peptides; the seven hollow markers are where a side chain was rebuilt. The figure shows that the binding site is shared — it is not independent evidence that BoPEP4 binds there.
[potentially add graph: the five AlphaFold seeds for free BoPEP4 as a per-residue pLDDT overlay, with the secondary structure each seed assigned shown as a strip beneath, purpose: to make fraction_disordered = 1.0 visible as five disagreeing predictions rather than a number, which is what motivated threading]
3.2 Receptor and co-receptor structures
Two receptors and one co-receptor make up the BoPEP4 signalling complex. PEPR1, an LRR receptor-like kinase acting as an innate immunity pattern recognition receptor, has its extracellular domain in 5GR8 chain A, which we used as the docking receptor input.
Given the similarity between PEPR1 and PEPR2 we considered PEPR2 as an alternative. PEPR2 is more specific toward AtPEP ligands where PEPR1 shows broader recognition. We generated a PEPR2 structure with both AlphaFold 3 and SWISS-MODEL homology modelling against PEPR1, but selected PEPR1 as primary on the strength of its experimental structure, using PEPR2 as a complementary model.
BAK1 came from 4MN8 chain B, the FLS2–BAK1 complex bound to flg22. BAK1 is a shared LRR co-receptor for multiple plant immune receptors including the PEPRs, and its conserved ectodomain provides the framework for receptor–co-receptor association. 4MN8 does not contain PEPR1, but the resolved BAK1 ectodomain transfers on that conserved architecture.
[Figure: the three receptor-side molecules side by side, annotated with provenance — PEPR1 solved, PEPR2 modelled, BAK1 borrowed from another complex]
3.3 Receptor restraint mapping
Restraints mark the residues expected to participate in an interaction, constraining the search toward interfaces that are physically and biologically plausible. Given the limited and heterogeneous evidence for the molecules in this complex, we cross-validated across several tools before writing anything into a restraint file.
MUSCLE alignment and WHISCY
To supply WHISCY with evolutionary information for PEPR1 we collected homologues by BLAST and aligned 352 sequences with MUSCLE. WHISCY weights each surface residue by conservation, on the reasoning that conserved solvent-exposed patches are more likely to be functional interfaces. We ran it both with and without the alignment so the effect of the evolutionary input could be seen rather than assumed.
Results · what the scan and the predictor returned
41 active residues, 7 supported by both
The primary restraint set is contact-derived from 5GR8: every PEPR1 heavy atom against every AtPEP1 heavy atom at 4.5 Å, across both crystallographic copies, unioned. That gave 41 active residues covering most of the concave LRR face. Passive residues were solvent-exposed positions (relative SASA ≥ 0.40) with a Cα within 8 Å of an active residue, giving 9. The passive set is thin because the active interface is unusually extensive, as expected for an LRR wrapping its ligand.
The same scan found the anchor Tang et al. report, without being pointed at it: Arg487 salt-bridges the C-terminal carboxylate of Asn23 at 2.40 Šand 2.59 Šin the two copies. Arg491 makes no hydrogen bond and no salt bridge on our PDBePISA run. It buries 12 Ų and sits about 5 Šaway, so we record it as an interface residue and not as a second anchor, and we do not test models against it.
On the peptide side the scan returned 16 of the 17 resolved AtPEP1 residues, and we treated the peptide as fully active, standard for a short linear ligand enclosed in a groove.
Run blind, WHISCY's active tier recovered little of the crystal interface, 1 of 41 and 2 of 41, agreeing with each other on a single residue, so we did not treat it as a reliable independent check. The passive tier behaved differently, returning the same six crystal-contact residues either way. Seven residues (347, 348, 369, 371, 372, 393, 395) appear in the crystal contacts and in both runs; we kept all 41 as active and flagged those seven as the highest-confidence core.
For PEPR2 we transferred the lists by pairwise alignment, giving 39 active and 9 passive. 21 of the 41 PEPR1 active residues are identical in PEPR2. Arg487, the anchor, is conserved as Arg464. These lists are a starting hypothesis for the homology model, not a validated restraint set.
[potentially add graph: the WHISCY and CPORT per-residue interface propensity scores for PEPR1 with the active/passive cut drawn, using the run output rather than a schematic, purpose: to show the restraint list was read off a threshold and not hand-picked, and to expose how many residues sat near the cut]
US-align and the BAK1 interface
No structure of PEPR1 bound to BAK1 exists, so the co-receptor interface had to be predicted by structural analogy with FLS2.
Results · the transferred interface
23 residues mapped onto PEPR1
FLS2 and PEPR1 share only 28% sequence identity, too low for sequence-based transfer, but US-align returns TM-score 0.840, well above the 0.5 threshold for the same fold. We read residue equivalences directly from the alignment block rather than approximating by nearest neighbour, then mapped FLS2's 25 BAK1-contact residues onto PEPR1, giving 23 predicted BAK1-site residues.
An initial placement rotating the whole BAK1 domain produced a physically impossible intertwined model. A local Kabsch refit restricted to matched interface residues brought interface RMSD from 4.00 Å to 2.09 Å and raised plausible contact pairs from 61 to 95.
Two of the 23, residues 372 and 395, fall inside the seven-residue peptide-interface core, and several sit near the 487–491 region. In the FLS2 system the ligand's C-terminal region and BAK1's footprint are structurally adjacent, so this is what the architecture would predict. We report it as consistent with known biology, not as evidence for it, since a geometric transfer could land generically nearby.
On the BAK1 side the restraints are contact-derived from 4MN8: 24 active and 4 passive, comprising contacts with FLS2 and with flg22, where Val54 and Thr58 are shared. We had expected Koller and Bent (2014) to give a second, conservation-based line of support. On reading it in full we found it reports the opposite: the conserved cluster it identifies did not impair FLS2–BAK1 interaction under mutagenesis, and the paper concludes FLS2 and BRI1 use receptor-specific docking sites. The BAK1 restraints therefore rest on a single line of evidence and were weighted accordingly.
[potentially add model: PEPR1 and PEPR2 superposed by US-align with the TM-score reported and the restraint-bearing residues marked on both, purpose: to justify transferring a restraint set derived on PEPR1 to PEPR2, which the page currently supports with an identity percentage alone]
DisVis
Before committing to full runs, DisVis tests whether a restraint set describes a geometrically possible arrangement at all. It scans every rigid-body placement of one partner around the other and reports how many restraints can be satisfied simultaneously, with the volume of accessible interaction space. Each line in a DisVis file is independent, with no OR grouping, so a broad active set is not expected to be fully satisfiable; what matters is how much of the set and how much space stay consistent, measured against a case where the answer is known.
Results · the two runs
Prepared, not yet run
Because we had active residue lists rather than known pairs, we generated the full combinatorial matrix for each run, as Cα–Cα distances bounded 2.0 to 8.0 Å. Run 1, PEPR1 against AtPEP1 (656 restraints), is the positive control drawn from a solved co-crystal. Run 2, PEPR1 against BAK1 (552 restraints), is the test: its two sides were derived independently and never checked for mutual compatibility.
DisVis has no numbers on it yet. The inputs are prepared, but the web server needs a WeNMR login we do not have, so the satisfiable-restraint count, the accessible interaction space, and the Run 2 against Run 1 comparison are all missing.
[potentially add graph: the DisVis accessible interaction volume as restraints are added one at a time, from the actual run rather than the illustrative version in the map, purpose: to show the final restraint set is satisfiable and by how much margin, which is the only quantitative evidence that the restraints do not contradict each other]
Contact-derived restraints from the crystal structures
Structural superposition gave reference geometry for each binary interface, which is what lets us judge docked models on structure rather than on score alone.
Results · the reference poses
fnat 0.974 on the threaded peptide
For the peptide interface we threaded BoPEP4 onto the AtPEP1 backbone from 5GR8 chain J. The gapless alignment means numbering carries over one to one. The threaded model reproduces the native interface with fnat 0.974, so it is the reference pose for BoPEP4 in the groove. That number is mostly a consequence of how the model was built: grafting onto the native backbone recovers the native contacts at conserved positions by definition. It tells us the threading came out undistorted. It does not tell us BoPEP4 binds this way.
Asn23 is held by five contacts and not by one salt bridge: OXT–Arg487 NH1 at 2.40 Å, O–Arg487 NE at 2.80 Å, ND2–Asp441 at 2.71 Å, N–Asp441 at 2.75 Å, and O–Asn465 at 3.14 Å.
For the co-receptor interface we superimposed FLS2 (4MN8 chain A) onto PEPR1 along the LRR β-ladder at 1.59 Å over 74 Cα pairs and carried BAK1 with the transform. As a check, the transformed flg22 lands within 1 to 2 Å of AtPEP1, confirming both ligands occupy the same register on their respective receptors.
HADDOCK score did not reliably track structural accuracy in these runs: one cluster ranked fourth by score returned fnat 0.000. We therefore rescored every cluster against the 5GR8 and 4MN8 reference geometry and selected models by agreement with those references rather than by score ranking.
[potentially add model: the 5GR8 interface with all 25 unique PISA interactions drawn as dashed bonds, coloured by hydrogen bond against salt bridge, purpose: to turn the interaction inventory from a count into something a reader can inspect, and to show where on the groove the polar network concentrates]
[Figure: the threaded BoPEP4 reference pose overlaid on AtPEP1 from 5GR8 chain J, with the C-terminal anchor contacts drawn and measured]
3.4 LEA14 structural input
We took the protein from 1XO8, the solution NMR structure of At1g01470, deposited by the Center for Eukaryotic Structural Genomics. Two features of that deposition matter. It is annotated UNKNOWN FUNCTION, so nothing in the record says where the protein binds anything. And the construct is the whole protein: the deposition cross-references UniProt O03983 across residues 1 to 151, and the coordinates carry all 151 with no gaps and no expression tag. What we dock is the complete mature protein.
Results · what the 20-model ensemble contains
A rigid core and one mobile tail
Over the first ten models the mean pairwise Cα RMSD across all 151 residues is 3.94 Å, maximum 5.55 Å. Restricted to the folded core, residues 24–145, the same comparison gives 0.66 Å, maximum 0.84 Å. The spread is almost entirely in the N-terminal tail; the β-sandwich core is one conformation reproduced twenty times.
Per-residue fluctuation agrees. Superposing all 20 models on the core, Cα RMSF averages 8.75 Å across residues 1–16 against 0.42 Å across 24–145. The transition is sharp rather than gradual: Leu16 fluctuates by 3.97 Å and Ala18 by 1.16 Å, a factor of three across two residues. The C-terminal six residues sit at 1.61 Å, mobile without being disordered. The most rigid position is Cys58 at 0.26 Å.
Those figures reproduce the ones already in the LEA14 design documents to the same decimal, which is the cross-check we wanted before building on them. They also set the truncation boundary: a cut at residue 16 removes the mobile segment and leaves the ordered region intact, where a cut at 23 would remove six residues belonging to the folded core. Glu23 is 10.9% solvent exposed on our own calculation, which is to say buried.
[Two truncation boundaries, and neither was built. The engineering axis cuts residues 1–23, from the published structure's 24–145; this ensemble puts the boundary at 16. We did not choose between them. Both went onto the September design list and into Level 0 assembly on 26 September, and neither reached a sequence-verified construct; every LEA14 plasmid we hold is full-length.]
For docking we prepared two forms: model 1 alone, and the first ten models concatenated as a HADDOCK ensemble input. Chain and segment identifier A, heavy atoms only.
3.5 AtPP2-B11 structural input
No experimental structure of AtPP2-B11 exists, so this one had to be predicted. We folded the full 257-residue sequence of UniProt Q949S5 with AlphaFold 3 on the AlphaFold Server, templates enabled, and the server returned five models from the submitted seed.
Results · confidence and what the five models mean
pTM 0.89, mean pLDDT 90.8
All five models report pTM 0.89 with ranking scores 0.89 to 0.90, no clashes, and a disordered fraction of 0.00 or 0.01. Mean pLDDT on the top model is 90.8: 92.8 over the F-box region, 89.9 over the linker, 90.6 over the PP2 domain. Three stretches fall below pLDDT 70: residues 1–2 at the terminus, and two loops inside the PP2 domain at 179–186 and 196–202.
As an independent check we compared the top model against the AlphaFold Database entry for the same protein, computed separately. The two agree to 1.27 Å Cα RMSD over all 257 residues, and both place the same two loops in their low-confidence regions.
The word ensemble needs qualifying. Mean pairwise Cα RMSD among the five models is 0.94 Å, maximum 1.14 Å, and restricting to the PP2 domain changes nothing at 0.95 Å. The five models are the same conformation five times: they sample the predictor's uncertainty, not the molecule's conformational space. This is the opposite of LEA14, where the deposited models carry genuine variety in the tail.
Which part we dock
InterPro places the F-box domain at residues 1–44 and the PP2 lectin-like domain at 92–255, with an unassigned linker between. The F-box is the surface that binds ASK1 inside the SCF complex, so it is not available to LEA14 in the cell. Docking the full chain would let the search place LEA14 on a face that is occupied in vivo, so the input is residues 45–257 with the F-box removed, numbering preserved. That leaves the PP2 lectin domain as the surface under test, which is also the substrate-recognition module by analogy with the solved plant F-box complexes.
Binary docking reduces the complexity of the full system and characterises one interface at a time. Most of our runs, including the alanine scans and truncation ladders, are built on this workflow because of its relative simplicity.
We docked BoPEP4 against the PEPR1 ectodomain in HADDOCK 2.4. The conformations generated during structure preparation were treated as alternative representations of the flexible peptide, and the active and passive residues from restraint mapping guided the search toward plausible interfaces.
Sampling was 1000 rigid-body, 200 semi-flexible and 200 water-refined models per run, FCC clustering at 0.60, minimum cluster size 4, peptide fully flexible. Confirmation runs, which ask the cost of a single edit rather than where the peptide docks, used 40/10/10 on the same protocol. Settings were constant across every binary run, so scores stay comparable.
Models were clustered by structural similarity and evaluated by HADDOCK score and subsequent structural assessment. We did not treat the highest-scoring structure as correct, and retained structurally plausible clusters for comparison.
Every cluster representative from the nine BoPEP4 runs, placed by its HADDOCK score and its interface RMSD from the crystal. Marker area is cluster population. A converged search produces a funnel that narrows towards the bottom left; these runs produce a cloud, and the two bands near 4–6 Å and 11–14 Å are the bound-like and inverted placements.
4.2 Threaded ternary docking
BoPEP4's activity is not determined solely by PEPR1. BAK1 is the signalling kinase balancing growth and innate immunity, and PEP perception involves receptor–co-receptor association, so we asked whether the predicted BoPEP4–PEPR1 interaction could be accommodated in a PEPR1–BoPEP4–BAK1 complex.
Rather than docking all three components simultaneously, we used the binary result as a structural starting point: a selected PEPR1–BoPEP4 complex was retained as established, and BAK1 was introduced on the interface identified during restraint mapping. This preserved the geometry of the modelled ligand–receptor interaction while searching for a compatible BAK1 orientation. Resulting models were assessed for steric clashes, interface plausibility and consistency with the restraints.
Two caveats. The ternary run inherits its binary input, so a peptide-register error would carry through uncaught, which is why the binary complex was validated against 5GR8 first. And its two restraint sets are not equally certain: peptide-to-receptor is contact-derived from a solved co-crystal, BAK1-to-receptor transferred from 4MN8 by analogy. The output is a structural hypothesis consistent with the published heterodimerisation biochemistry, not a residue-level prediction.
[potentially add model: the ternary PEPR1/BoPEP4/BAK1 assembly with the inferred co-receptor face shaded differently from the contacts taken from the crystal, purpose: to keep the inferred part of the model visually separate from the measured part, which the page insists on in prose and should enforce in the figure]
4.3 AI co-folding
HADDOCK searches under restraints we supply; a co-folding predictor builds from sequence with none. The two fail for different reasons, so agreement between them carries information. We ran AlphaFold 3 on the AlphaFold Server, templates enabled, five models per seed, with nothing from the HADDOCK workflow entering the job. The two tracks run in parallel and never seed each other.
We record ipTM, pTM, the ranking score and the chain-pair ipTM matrix with interface PAE, since a single ipTM averages across every interface and a confident pair can hide a failed one. Tracks are compared on interface residue identity, not global RMSD. Solved complexes of the same class are run first under identical settings, and the acceptance threshold is fixed from what they return.
Results · confidence is not activity
ipTM 0.73 for a peptide 4000-fold inactive
Our 30 July audit found AlphaFold 3 returned ipTM 0.73 for a peptide published as more than 4000-fold inactive, against 0.74 for BoPEP4. The metric does not separate a working ligand from a dead one, so we use it only as a structural confidence estimate. Multiple models from one seed are the predictor's uncertainty about one answer, not a conformational ensemble.
The docking and co-folding runs themselves are in Sprints and Results.
[potentially add graph: the predicted aligned error matrix and ipTM for each co-folding run, chains blocked out along both axes, purpose: to show inter-chain confidence directly, since ipTM alone does not reveal whether the uncertainty sits at the interface or inside a chain]
Docking returns a pose. It cannot say whether that pose survives once the atoms are allowed to move, and every design decision above rests on poses that have only been scored, not watched. Flex is where that gets tested.
The work is nine all-atom trajectories, including an AtPEP1 benchmark and an apo receptor control alongside the BoPEP4 constructs. It is reported in full on the Molecular Dynamics page: the systems and their parameters, what each run changed about the design, and the limits of a 30 ns window.
[potentially add graph: backbone RMSD against time for each replicate alongside per-contact occupancy bars, from the production trajectories, purpose: to show whether the complex equilibrated and which contacts persisted, which separates a real interface from a docking artefact]
Docking returns hundreds of models and will always return a best one, whether or not any of them are right. Triage decides whether a model is worth carrying forward, and it is deliberately built from methods that fail for unrelated reasons. A model passes only if it passes all five. When one fails, which one fails is the diagnosis that shapes the next round of truncations and substitutions.
Did the search converge?
Is the model actually right?
Is this the size and chemistry of a real interface?
Is the geometry chemically legal?
Does it reproduce contacts we already know exist?
6.1 HADDOCK score and clustering
HADDOCK groups similar models into clusters, each representing one binding mode. We read the funnel shape and the cluster population rather than the single best-scoring model, because a lone good score with no supporting cluster usually means the search sampled a lucky pose without converging. A well-populated cluster at the bottom of a funnel is the shape we want.
Why this is not fully applicable: the score discrepancy
HADDOCK score is a weighted sum of energy terms in arbitrary units, designed to rank models within a run. In our own validation it did not track structural accuracy: one cluster ranked fourth by score reproduced none of the native contacts. This is the single most important reason triage exists as a multi-check process. If the score alone were reliable we would not need the other four. We use it to judge convergence and to rank within a run, never as evidence that a model is correct.
It matters more again across constructs of different length, where a longer chain buries more surface and scores better for that alone. Down a truncation series we compare per unit of interface, never raw score.
The same 68 clusters, scored against the fraction of crystal contacts each one recovers. If the HADDOCK score ranked models by accuracy the cloud would run from upper left to lower right. Spearman ρ is −0.21 against fnat and +0.17 against L-RMSD: the correct sign in both cases, but far too weak to pick a model. This is the score discrepancy stated as a number.
6.2 CAPRI
CAPRI is the community standard for grading a predicted complex against a known answer, sorting models into acceptable, medium and high accuracy on contact recovery and ligand placement. It is the only check that measures correctness directly, and it only works where a reference structure exists. For us that is 5GR8. We grade every protocol change against it, so that any claim about a model with no reference rests on a protocol calibrated on one that has it.
One caution: the grade is an aggregate. A model can land in the medium tier with the single most conserved contact in the system far out of place. That is what the literature check below is there to catch.
All 68 cluster representatives against the two CAPRI criteria, with the acceptable and medium boxes drawn. Twenty-six clusters reach acceptable and three reach medium. None reaches high, and nothing comes close to the 1 Å L-RMSD that tier requires — the best L-RMSD anywhere in the set is 4.44 Å.
6.3 PDBePISA
PISA reports buried surface area, chemical composition, and whether a contact looks like a genuine assembly or a crystal packing artefact. It answers a question the docking score cannot: is this the right size and character for a complex that exists in solution? A model can score well while burying far too little surface to be real.
The 5GR8 reference list
We ran PISA on the crystal first, so every later number has a reference point. Our docked models are read against those values rather than a generic threshold, which matters for a 23-mer peptide where expectations built on globular interfaces do not transfer.
Results · the reference interface
2,039 Ų and 25 unique interactions
The 5GR8 interface buries 2,039 Ų, and that surface is 21% apolar, 41% polar and 38% charged. The bond inventory is 19 hydrogen bonds and 10 salt bridges, which comes to 25 unique interactions once the pairs that appear in both tables are counted once.
One property of this interface decides how we read PISA here: its own hydrophobicity test comes back negative. The peptide side gives a solvation gain of +2.5 kcal·mol−1 with a P-value of 0.585, so the interface is less hydrophobic than a random patch of the same size. The polar network is what holds it together, and that is why we read area and bond inventory here instead of the ΔiG verdict.
[potentially add graph: buried surface area and bond inventory for the crystal, the threaded reference and the best docked model as three grouped bars, purpose: to compare interfaces on the same axis, since the PISA numbers are currently reported for the crystal alone with nothing to compare them against]
6.4 MolProbity
MolProbity validates the stereochemistry of the coordinates themselves, independent of any docking score. It reports a clashscore for atoms in impossible overlapping positions, rotamer outliers for side chains in conformations real residues do not adopt, and backbone geometry. A model can satisfy every restraint and still be chemically illegal, particularly after the side-chain repacking refinement performs, which is why this check sits outside the energy-based scores.
The Ramachandran plot
The plot maps each residue's backbone dihedrals against the regions physically accessible to a polypeptide. Most residues should fall in the favoured regions, and outliers mark strained conformations. For threaded and refined peptide models this is a particularly direct check, because threading imposes a backbone taken from another sequence and this is where the strain from that imposition would show.
[potentially add graph: the Ramachandran plot of the refined complex with outliers labelled by residue, from the MolProbity run, purpose: to show the geometry check was passed by a specific structure, and to name any residue that failed it]
6.5 Literature check
The final check asks whether the model reproduces contacts that are already documented. For the PEPR1 interface that is Arg487 engaging the peptide's free C-terminal carboxylate, which we verified directly in the crystal coordinates before using it as a test. A model that fails to reproduce it is wrong regardless of how it scored on everything above, because it has placed the peptide in the wrong register.
We measured the criterion on the crystal before applying it to anything. Arg487 holds the carboxylate at 2.40 Å. Arg491 is often named alongside it, but on our PDBePISA run it makes no hydrogen bond and no salt bridge at all, so we report it as an interface residue and do not score against it.
This check is cheap, it is the hardest to fool, and it is the one that caught register errors the energy scores accepted.
[potentially add model: the C-terminal anchor in the crystal and in the refined model side by side, with all five clamp contacts drawn and their distances labelled on both, purpose: to show the literature check as a direct comparison, including the measurements that match the team's independent August numbers]
What triage cannot establish
A model passing all five checks shows the sequence is capable of the molecular recognition step. It does not show that the peptide reaches the receptor, survives the apoplast, or that binding produces the downstream immune output. Calcium influx, ROS burst, MAPK activation and receptor complex phosphotransfer are read out in the cell, not in the coordinates. The redesign loop optimises for predicted binding; only the wet lab closes the gap to function.
Which checks each sprint passed and which it failed is recorded run by run in Sprints and Results.
7Codon optimisation
7.1 ViennaRNA
We used ViennaRNA on the assumption that reducing mRNA secondary-structure stability around the MF001 ribosome binding site and the start codon improves translation initiation. Performance was assessed by minimum free energy and unpaired-probability predictions, with codon adaptation index values against Bacillus subtilis 168 codon usage. Recoding codons 2 to 31 within the SamyQ signal peptide and the early CDS breaks inhibitory structure in the start+15 and SD core regions, making the start site easier for ribosomes to attach to.
The approach follows Castillo-Hair et al. (2019) on optogenetic control of B. subtilis gene expression. All analyses ran on the RS-inclusive Window A, modelling the pre-recombination transient mRNA state: RS (28 nt) + MF001 RBS (24 nt) + spacer (1 nt) + SamyQ signal peptide (93 nt) + first 90 nt of CDS.
MFE
Minimum free energy (kcal/mol) of the most stable secondary structure.
start+15 / start+30
Mean unpaired probability over the first 15 or 30 nucleotides from the start codon.
SD core
Mean unpaired probability over the 7-nucleotide Shine-Dalgarno core (AAGGAGG).
CAI
Codon adaptation index relative to B. subtilis 168 usage from the Kazusa database.
Phase 1 held the 93-nt SamyQ sequence static and optimised the first 15 codons of each passenger gene independently across 8,000 random synonymous variants, maximising a balance score of start+15 plus SD core. Phase 2 addressed elongation via CAI and appended a C-terminal 6×His tag; this was the June and July design, and the tag moved to the N-terminus for BoPEP4 constructs ordered after 16 August. Full-length optimisation degraded initiation by lowering start+15. Phase 3 departed from standard host adaptation by making SamyQ codons 2 to 31 editable, further reducing structure around the initiation site.
[Figure 1 | Candidate variants (SamyQ + CDS editable), four panels: (A) BoPEP4, (B) LEA14, (C) ACCD, (D) AmilCP. The GR candidate in yellow is the Castillo-greedy candidate, prioritising start+15 at the expense of SD and CAI. The best candidate in green is the highest balance score.]
[Figure 2 | Predicted MFE secondary structures of RS-inclusive Window A before and after codon optimisation, four panels: (A) SamyQ-BoPEP4, (B) SamyQ-LEA14, (C) SamyQ-ACCD, (D) SamyQ-AmilCP. Each shows the 5′ end, the MF001 RBS, the start codon, and the mutated bases differing between original and optimised sequences.]
All four cargoes show the same pattern. The upstream RBS baseline barely moves, while the first stretch of the ORF rises in unpaired probability over a short window. Breaking the double-stranded stems at that point opens the start site, and the sequence in front of it never has to unwind.
[potentially add graph: per-position cleavage-site probability from SignalP for all four cargoes on one axis, with the predicted cleavage point marked on each, purpose: to show the signal peptide is cleaved where intended in every construct, not only in the one that is discussed]
7.2 SignalP 6.0
SignalP predicts whether a sequence carries a signal peptide, which secretion pathway it targets, and where signal peptidase cleaves. We use it for two purposes.
First, to confirm each signal peptide and mature sequence junction is cleaved where intended, so the released mature protein has the N-terminus we designed rather than a truncated or extended one. This matters because the mature N-terminal residues carry charge rules affecting secretion efficiency, and one of our substitutions sits in that region.
Second, to check that no cleavage site is cryptically predicted inside the mature sequence, which would cut the protectant itself. A construct that passes structural triage and fails this check is not orderable.
Predictions were run across the whole signal peptide panel rather than on a single assumed choice, because the best signal peptide differs for each cargo and cannot be predicted from first principles.
References
Tang, J., Han, Z., Sun, Y., Zhang, H., Gong, X. & Chai, J. Structural basis for recognition of an endogenous peptide by the plant receptor kinase PEPR1. Cell Research 25, 110–120 (2015). PDB 5GR8.
Pearce, G., Yamaguchi, Y., Munske, G. & Ryan, C. A. Structure-activity studies of AtPep1, a plant peptide signal involved in the innate immune response. Peptides 29, 2083–2089 (2008).
Sun, Y. et al. Structural basis for flg22-induced activation of the Arabidopsis FLS2 and BAK1 immune complex. Science 342, 624–628 (2013). PDB 4MN8.
Waese, J. et al. ePlant: visualizing and exploring multiple levels of data for hypothesis generation in plant biology. Plant Cell 29, 1806–1821 (2017).
Castillo-Hair, S. M., Baerman, E. A., Fujita, M., Igoshin, O. A. & Tabor, J. J. Optogenetic control of Bacillus subtilis gene expression. Nature Communications 10, 3099 (2019).
van Zundert, G. C. P. et al. The HADDOCK2.2 web server. Journal of Molecular Biology 428, 720–725 (2016).
Edgar, R. C. MUSCLE: multiple sequence alignment with high accuracy and high throughput. Nucleic Acids Research 32, 1792–1797 (2004).
de Vries, S. J. & Bonvin, A. M. J. J. CPORT: a consensus interface predictor. PLoS ONE 6, e17695 (2011).
Zhang, C., Shine, M., Pyle, A. M. & Zhang, Y. US-align: universal structure alignments. Nature Methods 19, 1109–1115 (2022).
van Zundert, G. C. P. et al. The DisVis and PowerFit web servers. Journal of Molecular Biology 429, 399–407 (2017).
Krissinel, E. & Henrick, K. Inference of macromolecular assemblies from crystalline state. Journal of Molecular Biology 372, 774–797 (2007). PDBePISA.
Williams, C. J. et al. MolProbity: more and better reference data for improved all-atom structure validation. Protein Science 27, 293–315 (2018).
Lorenz, R. et al. ViennaRNA Package 2.0. Algorithms for Molecular Biology 6, 26 (2011).
Teufel, F. et al. SignalP 6.0 predicts all five types of signal peptides using protein language models. Nature Biotechnology 40, 1023–1025 (2022).
Koller, T. & Bent, A. F. FLS2 and BRI1 receptor kinases use distinct docking sites for the co-receptor BAK1. Plant Physiology 164, 1196–1208 (2014).
Jumper, J. et al. Highly accurate protein structure prediction with AlphaFold. Nature 596, 583–589 (2021).
Abramson, J. et al. Accurate structure prediction of biomolecular interactions with AlphaFold 3. Nature 630, 493–500 (2024).
Sources
Tang J. et al. (2015) Structural basis for recognition of an endogenous peptide by the plant receptor kinase PEPR1. Cell Research 25:110–120. PDB 5GR8.
Yamaguchi Y., Huffaker A., Bryan A.C., Tax F.E., Ryan C.A. (2010) PEPR2 is a second receptor for the Pep1 and Pep2 peptides — source of the 76% PEPR1/PEPR2 identity figure.
Sun Y. et al. (2013) Structural basis for flg22-induced activation of the Arabidopsis FLS2–BAK1 immune complex. Science 342:624–628. PDB 4MN8.
Trellet M., Melquiond A.S.J., Bonvin A.M.J.J. (2013) A unified conformational selection and induced fit approach to protein–peptide docking. PLOS ONE 8:e58769.
Charitou V., van Keulen S.C., Bonvin A.M.J.J. (2022) Cyclisation and docking protocol for cyclic peptide–protein modeling. J. Chem. Theory Comput. 18:4027–4040.
de Vries S.J., Bonvin A.M.J.J. (2011) CPORT: a consensus interface predictor. PLOS ONE 6:e17695.
de Vries S.J., van Dijk A.D.J., Bonvin A.M.J.J. (2006) WHISCY: what information does surface conservation yield? Proteins 63:479–489.
van Zundert G.C.P., Bonvin A.M.J.J. (2015) DisVis: accessible interaction space of distance-restrained complexes. Bioinformatics 31:3222–3224.
Karaca E., Melquiond A.S.J., de Vries S.J., Kastritis P.L., Bonvin A.M.J.J. (2010) Building macromolecular assemblies by information-driven docking. Mol. Cell. Proteomics 9:1784–1794.
Krissinel E., Henrick K. (2007) Inference of macromolecular assemblies from crystalline state (PDBePISA). J. Mol. Biol. 372:774–797.
Koller T., Bent A.F. (2014) FLS2–BAK1 extracellular domain interaction sites required for defense signaling activation. PLOS ONE 9:e111185.
Yin R. et al. (2022/2023) Benchmarking AlphaFold-Multimer on multi-chain complexes; and the BRI1–BAK1 plant receptor-kinase test case, Bioinformatics (2023).
Zhang C., Shine M., Pyle A.M., Zhang Y. (2022) US-align: universal structure alignment. Nature Methods 19:1109–1115.
Edgar R.C. (2004) MUSCLE: multiple sequence alignment with high accuracy and high throughput. NAR 32:1792–1797.
Williams C.J. et al. (2018) MolProbity: improved all-atom structure validation. Protein Science 27:293–315; wwPDB validation pipeline.
Structural insights into the plant immune receptors PRRs and NLRs, Plant Physiology (2020) review — for the structure-to-function caveat.
Numerical values inside figures are illustrative of the expected shape of each result — funnel depth, percentile position, occupancy — not measurements from a completed run. Residue identifiers, PDB codes, protocol parameters and the BoPEP4 sequence are literal.
Before any modelling, decide what is worth modelling. Four stress axes are read in parallel and filtered down to peptides that are small, secretable, and have a receptor you can actually build.
The filters are ordered by how cheaply they eliminate a candidate. “Receptor identified” is last and hardest — it is also what makes the rest of this pipeline possible at all.
Data it uses
Stress-response literature per axis
Known protectant families
What it produces
A ranked candidate shortlist
One sequence promoted to Phase 02
Why this step exists
A docking pipeline will happily produce a beautiful model of a peptide that has no receptor, cannot leave the cell, or is never expressed under the stress you care about. The screen is cheap; the docking is not.
What each protectant was selected for
BoPEP4
Biotic stress, by priming the plant’s own alarm
A plant elicitor peptide: it is not itself protective, it tells the plant to defend itself. That makes its whole value conditional on a receptor existing and being reachable, which is why it is the one protectant that carries a docking pipeline at all.
LEA14
Osmotic and desiccation stress
A Late Embryogenesis Abundant protein, the family plants accumulate when cells dry out. It works by being present rather than by binding a receptor, so the screen asked whether it folds, whether it is stable, and what partner it associates with.
ACC deaminase
The stress hormone itself
ACC deaminase breaks down 1-aminocyclopropane-1-carboxylate, the immediate precursor of ethylene. Lowering ethylene reduces the growth inhibition plants suffer under drought and salt. It is an enzyme, so the screening question is catalysis and stability rather than binding.
Candidate set and stress axes as defined in the project methodology; receptor identification for the Pep family from Yamaguchi et al. and Tang et al.
Every protectant enters here as a one-letter sequence and nothing else. The job of this stage is to work out what kind of molecule it is, because that single answer decides how it gets coordinates, whether it can be docked at all, and which of the later stages apply to it.
Checks to run
Length — short enough for a groove to enclose, or a folded protein in its own right?
Net charge and whether the charged residues cluster
Hydropathy pattern — is there a membrane-sticky or aggregation-prone stretch?
Family membership, by alignment against relatives with known structures
Is there an experimental structure of this exact sequence, or of a close homologue?
What the answer changes
Short and disordered → thread onto the family's solved backbone rather than predicting one
Folds on its own → take the experimental structure if it exists, otherwise predict the monomer and validate it
Covalent chemistry → docking does not apply; inherit ligand positions by superposition
No receptor identified → the docking half of the map is unavailable, whatever else is true
What the three intakes returned
BoPEP4
A 23-mer that belongs to a family with a solved complex
Shaded tiles are basic; the dashed tile is the single acidic residue; the outlined tile is Asn23, whose free carboxylate is the anchor point in the receptor.
Same length class as AtPep1 and, like it, strongly cationic. The alignment is gapless with 10 of 17 ordered positions identical, which is what later licensed threading rather than prediction. The triage answer was: short, disordered, family member with a solved complex — so borrow the backbone and dock.
LEA14
A 151-residue protein that folds on its own
Not a ligand and not disordered. LEA14 is an intact domain with an experimental structure already deposited, 1XO8, so intake went straight to "take the solved coordinates" rather than to any prediction route.
The question for a protein that already has coordinates is which surface. A 446-sequence alignment scored every position for conservation and exposure, splitting the 151 residues into 18 conserved and exposed, 18 conserved and buried, and 115 neither. The conserved exposed set is where a partner would be expected to bind, and it resolves into two candidate patches plus a scatter of isolated positions.
Patch 1
Pro56, Ile57, Pro79, Asp125, Pro127, Gly130
Patch 2
Asp111, Trp112, Asp113
Triage answer: folded protein, structure available, no receptor in the PEPR sense — so the docking half of this map never applied to it.
ACC deaminase
A 338-residue enzyme, and a sequence check that caught an error
Our construct is Q5PWZ8, ACC deaminase from Pseudomonas putida UW4. The only crystal structure with both cofactor and substrate bound, 1TZ2, is a different organism (Pseudomonas sp. ACP, Q00740), so intake had to include a pairwise alignment between the two.
That check earned its place. Position 156 is tyrosine in the ACP sequence but phenylalanine in ours, so a substitution the design document had been calling Y156C does not exist in our construct. It is F156C. Three HADDOCK runs had already been submitted against the ACP background before this surfaced.
Intake also counted six cysteines with no disulfide partner — 41, 67, 95, 162, 196, 282 — which became the engineering target. Triage answer: folded enzyme, covalent cofactor chemistry, no docking.
Why this step exists
A pipeline will happily produce a confident model of a molecule it should never have modelled that way. Reading the sequence first is what separates the three routes: a disordered family member gets a borrowed backbone, a folded protein gets its own structure, and an enzyme whose ligands are covalently attached skips docking entirely. Getting this wrong is not caught later — every downstream number simply inherits the error.
Tang et al. 2015, Cell Research (extended AtPep1 binding mode); PDB 1XO8 (At1g01470); UniProt Q5PWZ8 and Q00740; Edgar 2004, NAR (MUSCLE).
Two design axes feed this stage: shorten the molecule, or change a residue in it. Both produce new sequences that re-enter at intake and earn their own run. Neither produces a result until something downstream survives the checks.
Truncation
Removes a liability — a hydrophobic stretch, a secretion problem, synthesis cost
Also shrinks the search space, which makes every later stage cheaper
Constraint: find the load-bearing end first and trim from the other one
Substitution
Adds or removes a specific interaction, or removes a chemical liability
Needs a reason that survives the structure, not just a plausible story
Constraint: a conserved position is conserved for a reason that may not be visible in your model
What each protectant actually proposed
BoPEP4
Truncation survived; substitution did not
The C-terminus is off limits: Asn23's free carboxylate is the anchor, and deleting the last residue is documented to compromise binding. So truncation went to the N-terminal side, where the hydrophobic GILIG block is a secretion liability and contributes nothing to the anchor. Three constructs were designed and ordered: 7–23, 9–23 and 15–23.
The substitution axis was opened and then closed. Two point mutations were carried for weeks before being withdrawn:
K18R was designed to create a salt bridge and a π-stack that an arginine makes in the crystal. Measurement showed wild-type Lys18 already reaches Asp348 at 2.62 Å and sits in cation–π range of Phe371. The gap it was built to fill was not there.
K8Q was never a binding design at all — it is a secretion lever, aimed at the positive-inside rule. It was retired from the binding axis as unmeasurable and deprioritised on its own axis.
A 2008 alanine scan of AtPep1 settles the general case: the positions where BoPEP4 differs from AtPep1 are alanine-tolerant, and the three that matter — Ser15, Gly17 and the free C-terminal Asn23 — are conserved. There is no structural case for improving BoPEP4 by point mutation. Truncation is the axis that survived.
LEA14
Surface patches identified, design frozen before it closed
With no receptor to bind, mutation ideation for LEA14 ran off conservation and exposure rather than off an interface. The 446-sequence alignment produced 18 conserved exposed residues, which cluster into two candidate patches — Pro56/Ile57/Pro79/Asp125/Pro127/Gly130, and Asp111/Trp112/Asp113 — with nine isolated positions and a flexible N-terminal tail at 1–16.
Those patches are hypotheses about where a partner binds, not validated sites. The strand was cut short before the design list was frozen, so no LEA14 variant reached an order.
ACC deaminase
Removing a liability rather than adding an interaction
The target is formulation, not catalysis. Six free cysteines oxidise and cross-link in a protein held in solution, which costs yield in the reactor and activity in the bottle. The design axis is therefore cysteine removal at six positions, each tested as alanine and as serine, with the family's own residue added at two of them.
Two substitutions were inherited from an earlier design document and carried through unchanged so the screen would test them rather than assume them: F156C and V293T. Both failed. Conservation put Phe156 at grade 8 with 97% aromatic across the family, and Val293 at grade 9 with valine at 99%; the stability screen then ranked V293T the worst of all seventeen variants at −2.71 kcal·mol⁻¹.
What survived: C41V, C67A, C162A, C282A, with Cys95 and Cys196 retained because the family consensus keeps a cysteine at 95 and nothing passes the cutoff at 196.
Read this as a hypothesis generator
These are design hypotheses, not predictions. Each variant is a fresh sequence: it goes back to intake, gets its own coordinates, its own run and its own score. A variant is only better if it beats the parent on the same pipeline — and ideally in an assay. Two of the three strands here ended by removing variants from their list rather than adding them, which is the more common outcome and the more useful one to record.
Tang et al. 2015 (C-terminal deletion sensitivity); Pearce et al. 2008 (AtPep1 alanine scan); Murakami et al. 1997 (C162A); ConSurf run 1791277560; DDMut job 17915452271577165.
The three protectants arrived with very different amounts of structural evidence, so each needed its own route to coordinates. One already had an experimental structure, one had a confident prediction, and one had no usable conformation at all.
The same question asked three times with three different answers. Only LEA14 has an experimental structure of itself; ACC deaminase has a high-confidence prediction; BoPEP4 has neither, which is why it is the one that had to borrow a backbone.
How each one got coordinates
BoPEP4
Two routes attempted, one kept
AlphaFold web server on the free peptide returned fraction_disordered = 1.0, pTM 0.18–0.35, and five seeds giving five different secondary structures
That is the expected result for a PEP-family peptide, which is disordered in solution and only orders on binding — so the free-peptide route gave nothing to carry forward
Threading onto AtPep1 from 5GR8 chain J: gapless alignment, 10 of 17 ordered positions identical, backbone transferred unchanged and side chains rebuilt at the seven differing positions
Why threading rather than an ensemble
The LRR groove is rigid and preformed while BoPEP4 is fully disordered, so the bound conformation is a property of the receptor more than of the peptide. Inheriting a backbone from a peptide already solved in that groove reproduces the moulding that actually happens; sampling three free-peptide geometries does not.
The transfer is also checkable, which sampling is not. The prolines BoPEP4 introduces at positions 10 and 13 sit at φ −81° and −85° on the borrowed backbone, both within proline’s range — so this sequence can physically occupy the backbone it inherited.
LEA14
The only experimental structure in the project
1XO8, the NMR solution structure of At1g01470, deposited as a 20-model ensemble. Nothing had to be predicted and nothing had to be borrowed.
LEA14 was also co-folded with its partner in AlphaFold 3 across five seeds, in two setups. That strand was cut short before docking, so the structures exist and the complexes do not.
ACC deaminase
A high-confidence prediction of our own sequence
The AlphaFold model of Q5PWZ8, 338 residues, mean pLDDT 97.7, exported through ConSurf so that conservation grades and coordinates share one numbering. It is the first structure of our own enzyme the project had.
Worth separating two things: the model is of our construct, but the ligands in it are not. PLP and ACC were transferred from a crystal of a different organism by superposition, so they are inferred and the file says so.
Tang et al. 2015, Cell Research (5GR8, extended binding mode); PDB 1XO8, solution structure of At1g01470; ConSurf run 1791277560 (ACCD model export).
The one piece of the assembly that needs no inference: a crystal structure of the PEPR1 LRR ectodomain with AtPep1 bound in the groove, including the C-terminal anchor pocket.
Left: the solenoid in cross-section — the peptide (accent) runs along the concave inner face, N-terminus at the top of the groove, C-terminus at the anchor. Right: the anchor itself, where the peptide’s free carboxylate is clamped by a pair of arginines.
Preparation
Keep chain A (receptor LRR); note chain of the bound peptide for later validation
Strip waters, sulfates and other crystallisation additives
Decide on glycans — remove unless you model them explicitly
Check for gaps in the LRR before docking across them
What it gives you
Ground-truth receptor coordinates
A contact list at ≤ 4.5 Å to convert into restraints
A reference complex to grade your own docking against in Phase 05
The cross-validation you get for free
Computing contacts directly from 5GR8 reproduces the Arg487/Arg491–to–C-terminus interaction that the paper reports independently from mutagenesis. When your own contact extraction recovers a published contact you did not put in by hand, your pipeline is measuring the structure correctly.
Tang et al. 2015, Cell Research 25:110–120; RCSB PDB entry 5GR8.
No experimental structure of PEPR2 exists. It is built from PEPR1 — a defensible substitution, because the two share 76% amino-acid identity and bind the same peptide family.
Homology transfer is strongest exactly where you need it: the LRR core that forms the peptide groove. Confidence falls off at the termini and inter-repeat loops — regions you should not be placing restraints on.
Two routes
Classical homology — SWISS-MODEL or MODELLER with 5GR8 as template. Transparent, and the alignment is auditable.
AI prediction — AlphaFold2 / ColabFold on the ectodomain sequence alone, not the full-length receptor.
Quality gates
Per-residue confidence across the binding groove
Ramachandran and clash check before docking (Phase 05 tooling, run early)
US-align against 5GR8 — the groove must superpose
Why the isolated domain
Evaluation of AlphaFold-Multimer on a plant receptor-kinase test case found that submitting only the interacting domain gives higher-confidence models than submitting the full-length receptor, whose transmembrane and kinase regions add noise the predictor cannot resolve.
◇ Report this as a model. Every PEPR2 result downstream inherits this assumption. Where PEPR1 and PEPR2 results diverge, first ask whether the difference is biology or modelling error.
Yamaguchi et al. 2010 (PEPR2 as second Pep receptor; 76% identity); AlphaFold-Multimer evaluation on BRI1–BAK1, Bioinformatics 2023.
BAK1 is a shared co-receptor: the same small LRR ectodomain that completes FLS2, EFR and BRI1 also completes PEPR1. No PEPR–BAK1 structure has been solved, so the geometry is borrowed from the one that has.
In the solved complex, BAK1 does not contact the ligand’s buried core — it grips the peptide where it emerges from the receptor groove, acting as a molecular glue that completes the binding site. That contact position, not the exact residue identities, is what transfers.
What is solid
4MN8 gives real BAK1 ectodomain coordinates (chain B, LRR region)
The architecture — co-receptor completing the ligand pocket — is conserved across the SERK-dependent LRR-RKs
AtPep1 is documented to induce PEPR1–BAK1 heterodimerisation
What is inferred
That PEPR1 presents BAK1 the same face FLS2 does
The exact contact residues on the PEPR side
The relative orientation of the two ectodomains
The honest framing
The evidence for PEPR1–BAK1 association is biochemical, not structural. Your ternary model is a structural hypothesis consistent with that biochemistry — which is a useful thing to have, and a different claim from “this is the complex.”
Sun et al. 2013, Science 342:624–628 (4MN8); Tang et al. 2015, Cell Research (PEPR1–BAK1 heterodimerisation assay).
The one supporting tool every strand used. An alignment answers three different questions depending on what you point it at: where a residue sits in another sequence’s numbering, how much a position varies across relatives, and whether the sequence you are modelling is the sequence you think it is.
Schematic: the two receptors are 76% identical overall, and the conservation profile peaks across the ligand-binding face. The two arginines that clamp the peptide C-terminus occupy conserved columns — which is what licenses transferring the restraint.
Procedure
Align the PEPR1 and PEPR2 ectodomain sequences with MUSCLE
Map every PEPR1 contact residue onto its PEPR2 column
Record the PEPR2 residue number — do not reuse PEPR1 numbering
Flag any contact that lands on a gap or a non-conservative substitution
What it produces
A PEPR2 active-residue list in PEPR2 numbering
A list of positions where the two receptors genuinely differ — the interesting ones
Why this step exists
PEPR1 and PEPR2 are partially redundant in Arabidopsis, and comparing their binding competence for the same peptide is a real biological question. That comparison is only meaningful if both docking runs are restrained to the same structural positions.
Three different alignments, three different jobs
BoPEP4
Two alignments, both load-bearing
A 352-sequence MSA built with BLAST and MUSCLE supplied the evolutionary input for the blind interface predictor, on the reasoning that conserved solvent-exposed patches are more likely to be functional.
A second, much smaller alignment did the heavier lifting. Aligning BoPEP4 against AtPep1 returned the result that licensed threading: gapless, with 10 of 17 ordered positions identical, so every backbone atom has a direct counterpart and numbering transfers one to one. A third moved PEPR1 contact residues into PEPR2 numbering by column rather than by arithmetic, since a single indel shifts every downstream residue ID.
LEA14
446 orthologs, scored for conservation rather than transferred
No numbering needed moving, because nothing was being docked. The alignment exists to score each of the 151 positions: 446 UniRef50 orthologs aligned with MUSCLE 5, conservation reported as Henikoff-weighted Shannon entropy, and exposure measured separately by Shrake-Rupley against Tien 2013 maxima.
Crossing those two axes is what produces the result: 18 residues that are both conserved and exposed, against 18 conserved and buried and 115 that are neither. The conserved exposed set is where a partner would be expected to bind, and it resolves into two candidate patches.
ACC deaminase
The alignment that caught an error before the DNA order
A pairwise alignment of our construct Q5PWZ8 against the crystal structure’s organism Q00740 found that position 156 is tyrosine in the crystal and phenylalanine in ours. A substitution the design document had been calling Y156C does not exist in our protein; it is F156C. Three docking runs had already been submitted against the wrong background before this surfaced.
A separate 150-sequence alignment, built inside ConSurf with MAFFT, supplied the per-position conservation grades. Its internal control was that all five catalytic residues came back at grade 9, the maximum — which is how you know the alignment is biologically sensible before reading anything else from it.
One caveat travels with those grades. Glu295 scores only 37% glutamate and Leu322 only 37% leucine, and those two residues are what distinguish ACC deaminase from its closest relatives. The 150 homologues therefore span the wider fold-type II family rather than ACC deaminases alone, so conservation here measures structural importance across a fold, not function.
Edgar 2004, NAR (MUSCLE); Katoh & Standley 2013 (MAFFT, used inside ConSurf); Henikoff & Henikoff 1994 (sequence weighting); Tien et al. 2013 (reference SASA maxima); Yamaguchi et al. 2010 (PEPR1/PEPR2 identity).
Sequence alignment says two residues occupy the same column. Structural alignment says they occupy the same place. For a homology model you want both to agree before you trust a transferred restraint.
The superposition is the check that matters for the model: if the LRR groove of the PEPR2 model lands on the PEPR1 groove, transferred restraints point at the right surface.
What to superpose
PEPR2 model onto 5GR8 chain A — global fold check
Then the groove region alone — local check where it counts
BAK1 (4MN8 chain B) onto any SERK homologue you also model
What to record
TM-score and aligned length
RMSD over the groove residues specifically
Any transferred restraint whose Cα deviates markedly after superposition
Why global RMSD alone will mislead you
A long LRR solenoid can have an excellent global score while one loop of the groove is displaced — and that loop may be exactly where your restraint sits. Score the interface region separately from the whole chain.
Superposition, where each strand used it
BoPEP4
Transferring a co-receptor interface across 28% identity
FLS2 and PEPR1 are too divergent in sequence for a sequence-based transfer, but US-align returns TM-score 0.840, well above the 0.5 threshold for a shared fold. That is what justified reading BAK1-contact residues off FLS2 and mapping them onto PEPR1 by structural equivalence rather than by arithmetic.
ACC deaminase
Superposition in place of docking, behind a stated gate
This strand uses superposition for something the docking strand uses docking for: placing a ligand. The AlphaFold model of Q5PWZ8 was superposed onto 1TZ2 in PyMOL, and the threshold was written down before the run — catalytic Cα RMSD below 1.0 Å with all ten identities matching. It returned 0.36 Å and 10 of 10, so PLP and ACC could be copied across.
LEA14
Not reached
No superposition was needed. LEA14 has its own experimental structure, so nothing had to be borrowed from a relative.
Zhang et al. 2022, Nature Methods (US-align).
Contacts read off a crystal structure tell you where the peptide bound in that crystal. A conservation-based predictor gives you an independent second opinion — one that does not depend on the structure you are trying to reproduce.
Residues supported by both the structure and an independent predictor become active (their burial is enforced). Their solvent-accessible neighbours become passive — allowed at the interface but not required, which keeps the restraint set from over-determining the answer.
Active residues
Must end up buried at the interface
Violating one costs energy in the scoring function
Keep the list short and well-supported — over-restraining manufactures the answer you assumed
Passive residues
Surface neighbours of the active set
No penalty if they stay solvent-exposed
Give the docking room to find the real register
The BAK1 special case
With no PEPR–BAK1 structure, the BAK1-side restraints come from two convergent sources instead: the contact geometry transferred from 4MN8, and a repeat-conservation analysis that identified a conserved cluster of solvent-exposed residues near the LRR C-terminus shared by FLS2, EFR, PEPR1 and BRI1 — a generic SERK-docking signature. Two independent methods pointing at the same patch is the strongest evidence available here.
The same question without a crystal to answer it
BoPEP4
An independent predictor that mostly disagreed
Run blind on PEPR1, the active tier recovered 1 and 2 of the 41 crystal-contact residues in two runs that agreed with each other on only one position. The passive tier behaved better, returning the same six crystal-contact residues either way. Seven residues appear in the crystal contacts and in both runs, and those were flagged as the highest-confidence core without being used to narrow the restraint file.
LEA14
Conservation and exposure, which is this logic without a partner
With nothing to dock against, the surface question was asked directly: which positions are both conserved and solvent-exposed? That returns 18 residues, resolving into two candidate patches — Pro56/Ile57/Pro79/Asp125/Pro127/Gly130 and Asp111/Trp112/Asp113 — plus nine isolated positions. They are hypotheses about where a partner binds, and nothing has tested them.
ACC deaminase
Not applicable — replaced by conservation grading
There is no protein–protein interface to predict. The equivalent question, which positions are load-bearing, was answered by ConSurf instead: an evolutionary rate per position from a 150-sequence alignment, with the five catalytic residues all returning grade 9 as the internal control.
de Vries & Bonvin 2011, PLOS ONE (CPORT); de Vries et al. 2006, Proteins (WHISCY); Koller & Bent 2014, PLOS ONE (repeat-conservation mapping).
Before spending compute on thousands of docking models, ask a cheaper question: is this combination of restraints even geometrically possible? DisVis does a systematic six-dimensional scan and answers in minutes.
Each added restraint carves away positions the ligand cannot occupy. A healthy restraint set converges on a small region; an inconsistent one converges on nothing, and DisVis reports which restraint is the outlier.
What it catches
Two “contacts” on opposite faces of the receptor
A restraint transferred to the wrong residue number
A peptide too short to satisfy all restraints at once
What it will not catch
A geometrically satisfiable but biologically wrong interface
Restraints that are all consistent because they all came from the same wrong assumption
Where it earns its keep here
Most usefully on the BAK1 restraints, which are inferred rather than measured. If the conservation-derived patch and the 4MN8-transferred contacts cannot be satisfied simultaneously, you have learned something important before running the ternary docking, not after.
Where the pre-screen applies
BoPEP4
Two runs, one of them a calibration
The PEPR1–AtPep1 pair is drawn from a solved co-crystal, so it is the positive control: it shows what a true interface looks like under this restraint scheme. The PEPR1–BAK1 pair is the real test, because its two restraint sets were derived independently and had never been checked for mutual compatibility.
ACC deaminase
Not applicable — no bodies to place
DisVis asks whether a set of distance restraints can be satisfied by some rigid-body arrangement. With the cofactor and substrate covalently attached, there is no arrangement to search over and no restraint file to screen.
LEA14
Not reached
No restraints were ever written for LEA14, so there was nothing to pre-screen.
van Zundert & Bonvin 2015, Bioinformatics 31:3222–3224.
HADDOCK in its peptide configuration: the prepared protectant goes in as the flexible partner, the restraint list drives the search, and the sampling is scaled up because a peptide has far more to rearrange than a rigid protein.
The funnel narrows twice: once by energy after rigid-body sampling, once by clustering after refinement. The peptide protocol raises the flexible-stage sampling fourfold because a 23-residue chain has to reorganise far more than a folded partner does.
Run it twice
Peptide → PEPR1 (5GR8 chain A)
Peptide → PEPR2 (Phase 02 model, PEPR2 numbering)
Identical protocol, so the two scores are comparable
Sanity checks on the output
Is the top cluster populated, or a single lucky model?
Does the top cluster put the C-terminus in the anchor pocket?
Are restraint violations low in the top cluster?
Why compare the two receptors
PEPR1 and PEPR2 are partially redundant in planta. A peptide that docks well to both, or markedly better to one, is a testable prediction — provided you remember that the PEPR2 side carries the model’s uncertainty with it.
What docking produced, per protectant
BoPEP4
Nine runs, and none of them usable
Nine two-body runs were submitted across the wild type, point mutants and tagged constructs. Measured against the crystal, no cluster representative came within 4 Å of the native backbone and the best reached 4.44 Å. The run with the tightest convergence placed the conserved C-terminal anchor 13 Å out of position.
Two faults explain it. The nine runs each started from a different peptide conformation, so the series was never controlled. And every run built the termini neutral, deleting the C-terminal carboxylate that makes the 2.40 Å anchor contact. The strand moved to threading, and HADDOCK was kept for refinement only.
ACC deaminase
Not applicable, and this is the most important decision on the page
We know this tool, we had already run it three times on ACC deaminase, and running it per variant was the obvious next move. It would have produced numbers. The numbers would have meant nothing: docking searches for the best non-covalent pose, and in this enzyme both the cofactor and the substrate are covalently attached before chemistry starts.
The published structural work on ACC deaminase reflects that. Positions come from crystallography and superposition, not from docking. We followed the field rather than our own familiarity.
LEA14
Not reached
The strand was cut short before any docking run was set up.
Trellet et al. 2013, PLOS ONE; Charitou et al. 2022, JCTC; Yamaguchi et al. 2010.
Take the winning binary complex and dock BAK1 onto it as a third body. This is the step where the evidence is thinnest — and the step whose restraints therefore deserve the most scrutiny.
Two restraint sets with very different provenance hold this assembly together. Set A is measured; set B is inferred from a homologous complex and an independent conservation analysis. Reporting them as equally certain would be the mistake.
Procedure
Take the top-scoring peptide–receptor model from Phase 04
Add BAK1 (4MN8 chain B ectodomain) as body three
Apply BAK1 restraints from Phase 03
Keep the binary interface restrained so it is not disrupted
Read the result as
A structural hypothesis consistent with the heterodimerisation biochemistry
Not a prediction of the complex at residue-level accuracy
A source of testable contacts for mutagenesis
The one thing worth checking above all
In the solved FLS2 system, the co-receptor grips the C-terminal end of the ligand as it exits the groove. If your ternary model puts BAK1 somewhere else entirely — on the outer convex face, or against the N-terminal end — that is a red flag about the restraints, not a discovery.
Three bodies, three different meanings
BoPEP4
Receptor, peptide and co-receptor
The intended end point: the winning binary complex with BAK1 added as a third body, on an interface nobody has solved. It inherits the binary result, so a wrong register upstream propagates silently.
LEA14
A three-chain problem solved a different way
LEA14 also has a three-body assembly — itself, AtPP2-B11, and the SKP1-like adaptor ASK1 — but it was approached by AlphaFold 3 co-folding across five seeds rather than by multibody docking. Running the pair and the triple separately tests whether the LEA14 interface survives the adaptor being present.
ACC deaminase
A four-body assembly, built by superposition
ACC deaminase is a homotetramer, so it is the largest assembly in the project — and none of it was docked. The tetramer was built by superposing our model onto each of the four chains of 1TZ2. That assembly then answered a real question: Cys41 sits at the A–B interface and no other cysteine does.
Karaca et al. 2010, Mol. Cell. Proteomics (multibody docking); Sun et al. 2013, Science.
Run the same three chains through an AI co-folding predictor and compare interfaces. This is a second opinion from a method with completely different failure modes — not a replacement for the docking.
The PAE matrix shows the predictor’s confidence in the relative placement of every pair of chains. Off-diagonal blocks are the ones that matter: a pale receptor–co-receptor block means the model has no strong opinion about where BAK1 sits, which is exactly the region your docking is least certain about too.
How to compare
Superpose the AI model and the HADDOCK model on the receptor
Compare interface residue lists, not global RMSD
Score each chain pair separately with pDockQ2
How to read disagreement
Agreement → meaningful, independent support
Disagreement → a region needing scrutiny, not an automatic loss for either method
Both confident and both wrong is possible — especially for shallow peptide interfaces
Why neither method wins by default
Benchmarks of AlphaFold-Multimer, AF3 and ColabFold on protein–peptide complexes show accuracy varies sharply with interface type, and receptor–peptide–co-receptor assemblies are among the harder cases. HADDOCK, meanwhile, only knows what your restraints tell it. Their errors are uncorrelated, which is precisely what makes the comparison informative.
AlphaFold did something different for each
BoPEP4
A useful negative result
Run on the free peptide it returns fraction_disordered = 1.0, pTM 0.18–0.35, and five seeds giving five different secondary structures. That is not a failed prediction, it is the correct answer for a peptide that only orders on binding — and it is the reason the ensemble plan was abandoned for threading.
LEA14
The strand’s main computational output
Two co-folding setups, each across five seeds: LEA14 with AtPP2-B11, and the same pair with ASK1 added. This is where most of the LEA14 evidence sits, and it is independent of any docking because none was run.
ACC deaminase
The structure the whole strand is built on
No co-folding, but AlphaFold supplied the only structure of our own sequence that exists: 338 residues at mean pLDDT 97.7, exported through ConSurf so that grades and coordinates share one numbering. Every distance, burial figure and stability prediction downstream rests on it.
Docking produces a static pose. Molecular dynamics asks whether that pose survives being solvated, heated to physiological temperature, and left alone for tens or hundreds of nanoseconds.
Example trajectory shape. The anchor contacts persist for most of the simulation; a contact that survives only a fraction of the trajectory was a feature of the docking pose, not of the complex.
What MD adds
Solvent and ions the docking treated implicitly
Backbone relaxation the semi-flexible stage could not reach
A time-resolved measure of which contacts are real
What to report
Backbone RMSD and its plateau
Per-contact occupancy over the production run
Whether the peptide stays in the groove at all
Where it fits the ranking
MD is also how variants get separated when their docking scores are close: an affinity-based ranking (for example from umbrella sampling along a pulling coordinate) discriminates candidates that a single-pose score cannot. Expensive — reserve it for the shortlist that survives Phase 05.
Which strands have dynamics
BoPEP4
The only strand with trajectories
Used both to test whether a docked pose survives solvation and to separate variants whose static scores sit too close to call. A 30 ns run on the wild type ranked Lys8 nineteenth of twenty-three residues by contact occupancy, which is part of what retired the position-8 design.
ACC deaminase
Not performed, and the cost of that is worth stating
No simulation experience on the team and no time this season to acquire it, so every ACC deaminase structure is static. What that costs: no measure of flexibility, nothing about behaviour at the temperatures a sprayed droplet reaches, and no way to tell whether a substitution changes dynamics without changing stability.
The one dynamics-adjacent answer came free from the crystal. Comparing the four chains of 1TZ2 against each other shows which regions differ between copies — a weak proxy for mobility. It found the 39–50 lid loop is not mobile here, with mean deviation 0.11–0.15 Å against a whole-chain average of 0.23–0.25 Å.
LEA14
Not reached
No simulation was run. 1XO8 being a 20-model NMR ensemble does give some sense of which regions are well determined, which is the nearest thing available.
Standard all-atom MD practice; complementary to the HADDOCK water-refinement stage described in Trellet et al. 2013.
Rank by cluster, not by single model. The shape of the score-versus-RMSD plot tells you whether the docking found a real energetic minimum or merely a lucky pose.
A funnel — many models converging on both low RMSD and low score — is the signature of a solution the energy function actually prefers. A flat scatter with one good outlier is not.
Report per cluster
Cluster size (number of members)
Mean score and its standard deviation
Mean restraint-violation energy
Buried surface area
Warning signs
Top cluster with only 4–5 members
Score gap between clusters within the error bars
High restraint violation in the best-scoring cluster
Clustering cutoff matters here
The peptide protocol tightens the fraction-of-common-contacts cutoff to 5.0 Å. A looser cutoff merges genuinely different binding registers into one cluster and flatters the result; that is why this parameter travels with the protocol rather than being left at its protein–protein default.
What converged, and what replaced the score where it could not
BoPEP4
A score that did not track accuracy
The score is a weighted sum, 1.0 Evdw + 0.2 Eelec + 1.0 Edesolv + 0.1 EAIR, in arbitrary units. In our runs the best-scoring and most convergent cluster was also the one with the anchor 13 Å out of place, and its whole advantage sat in the electrostatic term.
A methodology audit later voided any claim resting on a score difference between two same-length peptides. For peptides of different length the prohibition is stronger still: a longer chain buries more surface and scores better for that reason alone, so the truncation series must be compared per unit of interface, not raw.
ACC deaminase
No funnel to read, so a different number with a stated error bar
Nothing was docked, so there is no cluster population and no score. Its place is taken by a folding stability prediction across seventeen substitutions, with a cutoff declared before submission and the method’s own blind-test error quoted alongside: 1.14 to 2.24 kcal·mol⁻¹, which is wider than most of the differences in the table. The ranking is usable; the individual values are not.
LEA14
Not reached
Nothing was scored.
Trellet et al. 2013, PLOS ONE; Charitou et al. 2022, JCTC.
You have a crystal structure of the PEPR1 binary complex. That makes this a rare case: you can grade your own prediction against ground truth using the community’s standard criteria, rather than trusting the score alone.
The tiers are nested: any high-quality model is also medium and acceptable. For a 23-residue peptide in a groove, reaching medium is a solid result and high is genuinely difficult.
References to grade against
5GR8 — ground truth for the peptide–PEPR1 binary complex
4MN8 — structural analogue for the co-receptor interface geometry only
What cannot be graded
The PEPR2 complex — no reference structure exists
The ternary complex — 4MN8 is an analogue, not the same assembly
Report these as ungraded rather than assuming the PEPR1 grade carries over
The useful calibration
Dock a peptide whose complex you already know — AtPep1 against PEPR1 — and see what tier your protocol reaches. That number is your pipeline’s accuracy on this system, and it is the honest error bar to attach to every novel sequence you run afterwards.
Only one strand can be graded at all
BoPEP4
Graded against 5GR8, and it is the grade that matters
This is the rare case where ground truth exists, so the protocol can be calibrated rather than trusted. The useful exercise is to dock a peptide whose complex is already known and see what tier the protocol reaches — that number is the honest error bar to attach to every novel sequence afterwards.
A caution the nine runs made concrete: a MEDIUM grade is an aggregate. One model reached fnat 0.538 while the single most conserved contact in the system was 13 Å out of place.
ACC deaminase
Not applicable, and should be recorded as such
There is no reference complex of this enzyme with a ligand docked into it, because the ligand is covalently bound and the positions were inherited by superposition rather than predicted. Quoting a CAPRI tier for a superposition would be a category error. The equivalent check is the superposition gate, which was stated in advance and passed.
LEA14
Not applicable — no reference and no prediction
No complex was predicted, and no experimental LEA14–PP2-B11 structure exists to grade one against even if it had been.
CAPRI assessment criteria as applied to peptide complexes in Charitou et al. 2022, JCTC.
What an interface analyser can tell you from coordinates alone: how much surface a contact buries, by what chemistry, and whether that is the size and character a complex existing in solution would have. It answers a question the docking score cannot.
How to use the numbers
Comparatively, between your own models run identically
As a sanity band — is the interface in the range real complexes occupy?
Together, because area catches implausibly small contacts that ΔG alone will not
How not to use them
As a stand-in for a measured affinity — buried area is not a Kd
Across models built by different pipelines
On a model whose top cluster was poorly populated
Ranked by raw score between molecules of different length
What each protectant got from it
BoPEP4
The reference interface, measured on the crystal first
PDBePISA on 5GR8 gives the numbers every model is read against: 2,039 Ų buried in total, composition 21% apolar / 41% polar / 38% charged, and an inventory of 19 hydrogen bonds and 10 salt bridges — 25 unique interactions once the pairs listed in both tables are counted once.
Two things that inventory settled. The anchor is a five-contact clamp on Asn23, not a single salt bridge, with the shortest contact in the whole interface being Asn23 OXT to Arg487 at 2.40 Å. And Arg491 makes no hydrogen bond and no salt bridge at all — it buries 12 Ų and nothing more, so it cannot be used as a pass/fail criterion the way the protocol had been treating it.
A caution specific to this system: PISA's own hydrophobicity test is negative here. The peptide side returns a solvation gain of +2.5 kcal·mol⁻¹ and a P-value of 0.585, meaning the interface is less hydrophobic than a random patch of that size. What holds it together is the polar network, so read area and bond inventory rather than the ΔiG verdict.
ACC deaminase
Used for a subunit question, not an affinity one
There is no docked complex to score, so the interface tool is pointed somewhere else entirely: ACC deaminase is a homotetramer, and a residue that looks exposed on one chain may be buried against a neighbour. Measuring burial on a monomer would then give the right number for the wrong molecule.
The finding changed the lead variant. Cys41 sits at the A–B subunit interface and no other cysteine does. Its 68 Ų exposure was a monomer artefact; in the tetramer it is 37 Ų, with the thiol itself at 15.4 Ų. Its partners are all in chain B: Phe13 at 3.85 Å in a sulfur–aromatic contact, Ser43 at 4.34 Å, Gly44 at 4.02 Å.
The two Cys41 side chains face each other across the two-fold axis at SG–SG 7.90 Å, far too distant for a disulfide — which independently confirms that all six cysteines are genuinely unpaired. Submitting 1TZ2 to PISA and checking whether Cys41 appears in the A–B residue list is the public cross-check on that.
LEA14
Not reached
The LEA14 strand stopped before any complex was scored. The co-folding runs with AtPP2-B11 exist and could be analysed this way, but no interface analysis was run and none should be reported.
A caution that applies to all three
Interface statistics were built mostly on protein–protein complexes with larger, better-packed contacts than any of these. A 23-mer in a groove buries less surface and leaves more of the site solvated; an enzyme's covalently bound cofactor is not an interface at all. Never compare raw buried area across molecules of different length — a longer chain buries more for that reason alone.
Krissinel & Henrick 2007, JMB (PDBePISA); PDBePISA run on 5GR8; tetramer measurements on 1TZ2.
A good score says two molecules fit together. It says nothing about whether the coordinates are chemically sensible. This check is deliberately independent of every energy function in the pipeline, and it is the same validation a depositor faces at the wwPDB.
What it reports
Clashscore — atoms overlapping by 0.4 Å or more, per 1000 atoms
Ramachandran — backbone angles outside the physically reachable regions
Rotamers — side chains in conformations real residues do not adopt
Hydrogens must be added first; clashscore is an all-atom metric
What does not apply
Rfree and real-space fit — a computed model has no experimental data
Saying so explicitly, rather than omitting those rows quietly, is part of reporting a model honestly
What it caught, per protectant
BoPEP4
Validating a borrowed backbone
Threading raises a specific risk: the backbone came from a different sequence, so any strain shows up here. Run against the crystal through an identical pipeline, the threaded complex scores a clashscore of 8.08 against 7.54 for 5GR8 itself — and 83 of its 90 clashes are internal to the receptor and identical to the crystal's. Threading introduced seven, all at rebuilt side chains.
Two of those were real model-building errors rather than physics. A proline built where the template had a glycine was placed on the wrong side of the nitrogen, giving a C(i−1)–N–CD angle of 40.6° against a normal 125°. And a histidine substituted for a valine collided with Ile319 at 0.74 Å overlap. Both resolved on refinement; neither would have been visible in a docking score.
The check also rules a mutation out. Gly17 sits at φ +82.5°, ψ −159.6°, a conformation essentially restricted to glycine, and the space a Cβ would occupy is already filled by Tyr276 and Val299. G17A is sterically blocked on two independent counts.
ACC deaminase
Validating a prediction before trusting anything built on it
Every distance, burial figure and stability prediction downstream rests on one AlphaFold model, so it is validated before use rather than after. Mean pLDDT is 97.7 across 338 residues, which is high enough to measure distances from, but a global confidence figure does not substitute for geometry: the per-position value at each of the eight residues being changed matters more than the mean.
The model then had to pass a stated gate before its ligand positions could be inherited from 1TZ2: backbone RMSD 0.61 Å over 331 residues, catalytic Cα RMSD 0.36 Å with a worst case of 0.58 Å, and 10 of 10 catalytic residues matching by identity. The threshold was written down before the superposition was run.
LEA14
The one case where validation is already done
1XO8 is an experimental NMR structure deposited as a 20-model ensemble, so it arrives with its own validation record and the spread between models is itself a measure of where the structure is well determined. No prediction needed validating, because nothing was predicted.
Why this is separate from the score
The energy function and the geometry validator disagree about what "good" means, and that is the point. A model can win on score while carrying a distorted backbone at the interface, and a distorted interface is exactly where a distorted score comes from. Running this check on the inputs as well as the outputs is cheaper than discovering the problem after the analysis.
Williams et al. 2018, Protein Science (MolProbity); wwPDB validation pipeline; ConSurf run 1791277560 (AlphaFold model export); PDB 1XO8.
The last gate and the cheapest one. Somebody has usually already measured something about your system. If the model contradicts a published measurement, the model is wrong, regardless of how it scored on everything above.
It is also the hardest check to fool. Every other test is an aggregate over a whole structure, and aggregates can be satisfied by the right total amount of the wrong thing. A single documented contact either is reproduced or it is not.
What the literature said, per protectant
BoPEP4
One hard contact, and one paper that closed a design axis
Measuring the crystal first is what separates a real criterion from an assumed one. Arg487 clamps the carboxylate at 2.40 Å; Arg491, which the protocol had treated as a second anchor, makes no bond at all.
The second literature result closed an entire design axis. Pearce 2008 ran an alanine scan of AtPep1(9–23) and found [A18] fully active in the alkalinization assay at 2.5 nM, ten times the EC₅₀. If deleting the whole side chain at position 18 costs under five to ten fold, then swapping lysine for arginine there must cost less — which retires K18R whatever the geometry says.
That paper was published in 2008 and sat unread until six weeks of docking had already been spent on the substitution axis. The lesson is cheap to state and expensive to learn: read the functional literature before modelling, not after.
ACC deaminase
A prediction that matched a 1997 measurement
One of the six cysteines has already been tested experimentally. Murakami and colleagues (1997) replaced Cys162 with alanine and found activity unaffected, while chemical modification of that same thiol abolishes activity completely.
Our own ranking placed C162A as the cheapest of the six Cys→Ala substitutions at −0.26 kcal·mol⁻¹, on conservation and stability grounds, before that paper was found. A computational ranking that agrees with an independent wet-lab result is the only validation in this strand that did not come out of a computer.
The same paper also shapes the target: a thiol reactive enough that modifying it kills the enzyme is exactly the liability worth removing, and removing it genetically costs nothing measurable.
LEA14
Not reached
No model was produced for LEA14 that could be checked against a published contact. The conservation analysis and the co-folding runs exist; the validation stage does not.
Do this one first, not last
Its position at the end of the list is about cost, not importance. Computing the documented contacts directly from the reference structure before docking anything is the cheapest way to confirm your measurement pipeline works at all — if your own extraction independently recovers a contact a paper reports, the pipeline is sound, and only then is it worth trusting on a structure that has no paper to check against.
Tang et al. 2015, Cell Research; Pearce et al. 2008 (AtPep1 alanine scan); Murakami et al. 1997, Biosci. Biotechnol. Biochem. 61:506.
Six molecules, three protectants, and four different kinds of evidence. This is where the honesty of the pipeline is set: what can be taken from a crystal of the right thing, what has to be built, and what is borrowed from a different complex or a different organism.
The bar under each molecule is how much of its final role is supported by direct experiment. Five of the six carry an inference of some kind, and the figure legend should say which.
Shared preparation
Strip waters, ions and crystallisation additives
Decide on glycans consistently across all of them
Check for chain breaks inside the binding region
Run geometry validation before use, not after
Charged N- and C-termini where a terminal carboxylate carries a contact
What it produces
Clean coordinate files on one standard
A contact list from 5GR8 at ≤ 4.5 Å for the BoPEP4 restraints
Ligand positions for ACC deaminase, transferred from 1TZ2 by superposition
A provenance note per molecule for the figure legend
What the partner side meant for each strand
BoPEP4
Three molecules, one of them solved
PEPR1 from 5GR8 is the only structure here of the exact molecule it stands for. PEPR2 is a homology model at 76% identity, and BAK1 is real coordinates taken from a complex with a different receptor entirely.
The preparation that mattered most turned out to be the least visible one. Nine docking runs were submitted with neutral termini, which omits the C-terminal carboxylate — the group that makes the shortest contact in the whole interface at 2.40 Å. No output model from any of the nine contained that atom, and nobody noticed for two months.
LEA14
A partner rather than a receptor
LEA14 does not bind a receptor, so its partner side is the F-box protein AtPP2-B11 (Q949S5, AlphaFold monomer, 257 residues), with ASK1 added in a second setup because an F-box protein does not function alone.
Neither has an experimental structure, and neither does the complex. Both the folds and the interface are predictions, which is the honest ceiling on anything this strand could have concluded.
ACC deaminase
A template, not a docking input
1TZ2 is here for a different reason from the others. ACC deaminase is a PLP enzyme: the cofactor is covalently bound to Lys51 as an internal aldimine and the substrate displaces it to form an external aldimine, so there is no pose to search for.
1TZ2 supplies PLP and ACC positions by superposition onto our own model, behind a gate stated in advance. That gate passed: backbone RMSD 0.61 Å over 331 residues, catalytic Cα RMSD 0.36 Å, and all ten catalytic residues matching.
◇ Report the transferred ligands as inferred. They are as good as the superposition and no better.
Tang et al. 2015 (5GR8); Sun et al. 2013 (4MN8); Yamaguchi et al. 2010 (PEPR2 identity); Karthikeyan et al. 2004, Biochemistry (1TZ2); UniProt Q949S5 (AtPP2-B11).
LEA14 was not modelled against a receptor. It was co-folded with AtPP2-B11, an F-box protein of the phloem protein 2 family, UniProt Q949S5.
What we have
AlphaFold monomer model, 257 residues
Co-folded with LEA14 (151 residues) in AlphaFold 3
Five seeds per setup, so the spread between seeds is readable
What is inferred
No experimental structure of AtPP2-B11 exists
No experimental structure of the LEA14–PP2-B11 complex exists either
Both the monomer fold and the interface are predictions
Where this strand stopped
The LEA14 line was cut short before it reached docking, so this molecule appears in the build stage and in the AI co-folding lane and nowhere else on this map. The co-folding runs exist and can be reported; the docking, restraint mapping and scoring stages were never run for LEA14.
UniProt Q949S5, F-box protein PP2-B11; AlphaFold 3 jobs S3_LEA14_AtPP2B11 and S3b_LEA14_AtPP2B11_ASK1.
An F-box protein does not function alone. The second LEA14 setup adds ASK1, the Arabidopsis SKP1-like adaptor that F-box proteins bind through, making the co-fold a three-chain problem rather than a pair.
Why it was added
F-box proteins work as part of an SCF complex
A two-chain prediction may place the F-box domain in a conformation the adaptor would not allow
Running both setups shows whether the LEA14 interface survives the addition
How to read the pair of runs
Interface unchanged between the two → the LEA14 contact does not depend on the adaptor
Interface shifts → the two-chain run was the artefact
Compare per-chain-pair confidence, not global scores
AlphaFold 3 job S3b_LEA14_AtPP2B11_ASK1, five seeds.
The only crystal structure of ACC deaminase with both cofactor and substrate bound. It is not our enzyme: 1TZ2 is from Pseudomonas sp. ACP, UniProt Q00740, while our construct is Q5PWZ8 from Pseudomonas putida UW4.
What it supplies
PLP and ACC positions, transferred to our model by superposition
Four chains in the deposited file, each with its own PLP and ACC
The tetramer geometry used for the subunit-interface check
The gate it had to pass
Backbone RMSD 0.61 Å over 331 residues
Catalytic Cα RMSD 0.36 Å, worst case 0.58 Å
10 of 10 catalytic residues matching by identity
The organism difference is not cosmetic
Checking it caught an error that had nearly reached the DNA order. Position 156 is tyrosine in the ACP sequence but phenylalanine in Q5PWZ8, so the substitution the design document called Y156C does not exist in our construct; it is F156C. Three earlier HADDOCK runs were submitted against the ACP background before this was found.
Karthikeyan et al. 2004, Biochemistry 43:13328 (1TZ2); UniProt Q5PWZ8 and Q00740.
The hub of the whole pipeline. Four supporting tools, none of which docks anything, converge here into a single artefact: the list of residues the docking is allowed to bury — and the confidence you have in each one.
The four supporting tools do not run in sequence — they run in whatever order is convenient and deposit evidence into the same table. The consistency scan is the gate that decides whether the table is worth docking with.
What the four tools contribute
MUSCLE — moves a residue identity from PEPR1 numbering to PEPR2 numbering
US-align — confirms the transferred position is structurally where it should be
WHISCY / CPORT — an opinion that does not come from the crystal
DisVis — does not add restraints; it rejects impossible combinations
What it produces
Active residues per molecule, with a source tag each
An ambiguous-restraint file the docking engine reads directly
Why this is the highest-leverage module
Information-driven docking cannot find an interface you did not point it at. Every downstream number — score, ΔG, CAPRI tier — is conditional on this file being right. It is also the cheapest thing to get wrong silently: a restraint transferred to the wrong residue number produces a confident, well-scored, completely fictional complex.
Who needed a restraint file
BoPEP4
41 active residues, 9 passive, 7 of them triple-supported
Contact-derived from 5GR8 at a 4.5 Å cutoff across both crystallographic copies. The passive set is thin because the active interface is unusually extensive, which is what an LRR domain wrapping its ligand looks like. Seven residues appear in the crystal contacts and in both blind predictor runs and were flagged as the core.
One property of these restraints explains a lot of what went wrong. An ambiguous restraint asks a peptide residue to approach any member of a 40-residue set. That constrains proximity and leaves register free, so several registers satisfy the file equally and the energy function picks among them.
ACC deaminase
Not applicable
Nothing is docked, so no restraint file exists. The stage that replaces it is structure preparation: building a model of our own sequence, validating it, and transferring ligand positions from a crystal behind a pass/fail gate.
LEA14
Not reached
No restraints were derived.
de Vries & Bonvin 2011 (CPORT); de Vries et al. 2006 (WHISCY); van Zundert & Bonvin 2015 (DisVis); Koller & Bent 2014 (conservation mapping).
One engine, two modes. Binary and ternary are not consecutive steps in a list — they are two runs of the same protocol with different numbers of bodies, and only a handoff of the winning model connects them.
The two docking modes share an engine, inputs and protocol; what differs is how many bodies are on the table. The AI lane consumes the same three chains but contributes nothing to the docking — keeping it independent is what makes agreement between them meaningful.
Shared by both modes
The prepared protectant as the flexible partner
The restraint table from Phase 03
Peptide-tuned MD steps and the 5.0 Å clustering cutoff
The same three-stage funnel: rigid body, semi-flexible, water
What differs
Binary: two bodies, restraints entirely from measured contacts
Ternary: three bodies, one interface inferred rather than measured
Ternary inherits the binary result, so a bad binary model poisons it silently
The dependency worth being explicit about
Ternary docking starts from the top-scoring binary complex. That is a single point of failure: if the peptide sits in the groove in the wrong register, BAK1 will be docked confidently onto a wrong assembly and nothing downstream will notice. Validate the binary complex against 5GR8 before spending a ternary run on it.
One engine, and two strands that never used it
BoPEP4
Both modes, plus a third use nobody planned
Binary and ternary ran as designed. The engine then found a third role it was better suited to: refinement only, with sampling switched off entirely, used to relax a threaded complex in explicit solvent without ever re-placing the peptide. That is the configuration that finally produced usable models.
ACC deaminase
The stage is deliberately empty
Three HADDOCK runs on ACC deaminase exist from before this was settled. They are kept in the record and not used, because the chemistry they assume is the wrong chemistry. An empty stage with a reason attached is a result.
LEA14
Not reached
The strand ended at structure generation and co-folding.
Trellet et al. 2013; Charitou et al. 2022; Karaca et al. 2010 (multibody).
Five checks that fail for unrelated reasons, run in any order, resolving to one verdict. It is a conjunction rather than an average — and which check fails is the diagnosis that writes the next design.
Each catches something different
Score & clustering — did the search converge at all?
CAPRI — is it right, measured against a real structure?
Affinity & interface — is this the size and chemistry of a real contact?
Geometry — are the coordinates chemically legal?
Literature — does it reproduce what is already known?
What the verdict feeds
A pass → report, with the provenance flags carried through
A fail → a specific diagnosis, not just a rejection
Either way → the next round of truncations and substitutions
How many of the five each protectant could actually run
BoPEP4 — all five
And the fifth overturned the other four
This is the strand where the conjunction earned itself. A set of nine docking runs produced a variant whose top cluster held 169 of 200 models, carried the best score in the series, and graded CAPRI MEDIUM at fnat 0.538 — the best grade anything in the set reached. It passed checks one, two and three.
Check five then put its C-terminal anchor 13.0 Å from Arg487, in a hydrophobic subsite lined by Leu392, Leu393, Ile415 and Thr417. The convergence was convergence onto the wrong register. Had the score alone decided, that model would have been the lead.
A second, blunter failure sits upstream of all five: across 68 cluster representatives, none came within 4 Å of the crystal backbone. The poses were wrong before any of them were scored.
ACC deaminase — three of five
Two checks do not apply, and saying so is part of the method
There is no docking run, so there is no funnel to read and no cluster to populate: check one does not apply. There is no reference complex of this enzyme with a ligand docked into it, so check two does not apply either — and claiming a CAPRI grade for a superposition would be a category error.
What remains runs in full. Interface analysis answered a subunit question rather than an affinity one. Geometry validated the model and gated the ligand transfer. The literature check matched a 1997 measurement.
In their place the strand added checks the docking pipeline does not have: an evolutionary conservation grade per position from a 150-sequence alignment, and a folding stability prediction with a stated cutoff and a stated error bar. The internal control on the first was that all five catalytic residues came back at grade 9 — the maximum — which is how you know the alignment is biologically sensible before reading anything else from it.
LEA14 — none
Stopped before scoring
The LEA14 strand reached structures and conservation analysis and went no further. Nothing was docked, so nothing was scored. The honest report is the two things that do exist: an experimental structure, and a conservation and exposure map with two candidate surface patches.
Read the failure, not just the flag
The point of five independent checks is that which one fails tells you what to change. A bad cluster means the restraints or the sampling. A good score with a missing documented contact means the wrong binding register. Clean geometry with a tiny interface means the sequence probably does not bind. Each sends a different variant back to the start — and a check that does not apply to your molecule should be recorded as not applicable rather than quietly skipped.
Charitou et al. 2022 (CAPRI on peptide complexes); Williams et al. 2018 (MolProbity); Tang et al. 2015; ConSurf run 1791277560; DDMut job 17915452271577165.