Journal · Simulation
Replicating a published protein folding result
One laboratory value matched, one of the toolkit's own values matched provisionally, and the test of whether the method generalises failed.
Updated 30 October 2026. The 1SHF result now comes from a run six times longer, the 3CLA comparison is marked provisional and identified as the toolkit's own simulation output rather than a laboratory value, and we say plainly below what the 1SHF test could and could not show. An independent review found these; the corrections are ours.
We took a published molecular simulation toolkit we had never run, on a laptop, and tried to reproduce numbers that laboratories, and the toolkit's own authors, had already published. This is the short version. The full record runs to about 153 minutes; this runs to about six.
What we compared, and how it came out
| Protein | Published | Ours | z | Verdict |
|---|---|---|---|---|
| 1SHF, Fyn SH3 domain | −6.0 ± 0.6 kcal/mol, laboratory | −6.039 ± 0.069 | −0.06 | agrees, but see below on what that shows |
| 3CLA, melting point | 341.2 K, toolkit output | 344.03 K | n/a | 2.83 K apart, no published uncertainty |
| 3CLA, stability | −8.4315 kcal/mol, toolkit output | −7.415 ± 0.859 | 1.18 | agrees, provisional: thin sampling, run never mixed |
| CI2, stability | −7.03 ± 0.16 kcal/mol, laboratory | −2.760 ± 0.021 | 26.5 | fails |
One laboratory value matched and one toolkit value matched provisionally. One failed, and the one that failed is the one we designed to test whether any of this generalises.
The result we are pleased with
The first row is not a fit. It is, though, a weaker test than we first said.
These models have one tunable knob, the strength of the forces holding the folded shape together. Tune it until the answer matches and you have described a measurement, not predicted one. So we measured the relationship instead: four separate runs at four different strengths, each one a full replica-exchange calculation.
Three of those points looked like a straight line. The fourth showed the line was wrong: the residuals of a straight line come out +0.182, −0.186, −0.183, +0.187, which is a curve, not noise. The curve as it stood on 23 October said a strength of 1.882 should land on the published value. Correcting an error in one calibration point two days later moved that to 1.890.
We then built a force field at 1.882, which had never existed before, and ran it. It returned −5.979 ± 0.115 kcal/mol against a laboratory value of −6.0 ± 0.6, a z of 0.03. Run six times longer it gave −6.039 ± 0.069, a z of −0.06. 1.882 was not one of the four points the curve was fitted to.
Here is what that does not show. The two calibration runs either side of 1.882, at 1.7453 and 2.0322, already agree with −6.0 on their own, at −5.303 and −6.597, because the laboratory uncertainty of 0.6 is wider than the whole gap between them. Any result between those two points would have agreed. So this confirms that the curve is smooth where we said it was. It was never a test the model could fail. The one that could fail is next.
The result that matters more
A force field that predicts the protein it was calibrated on is a much weaker claim than it sounds. The real question is whether the calibration carries to a molecule it has never seen.
So we took the strength calibrated on the SH3 domain, 1.882, applied it unchanged to chymotrypsin inhibitor 2, a different protein with its own published stability, and ran it. It does not transfer.
| Folded-state threshold | Our stability | z against published |
|---|---|---|
| Q ≥ 0.65 | −4.399 ± 0.037 | 16.0 |
| Q ≥ 0.70 | −3.596 ± 0.030 | 21.1 |
| Q ≥ 0.7321 | −3.018 ± 0.026 | 24.8 |
| Q ≥ 0.75 | −2.760 ± 0.021 | 26.5 |
| Q ≥ 0.80 | −1.959 ± 0.014 | 31.6 |
| Q ≥ 0.85 | −1.181 ± 0.010 | 36.5 |
The published value is −7.03 ± 0.16 kcal/mol, read from the primary source rather than a database. At the transferred strength the model gives CI2 about 39 per cent of its real stability, and the disagreement is far outside any uncertainty we can construct.
Why six rows rather than one. To compute stability you must decide which configurations count as folded, and earlier in this project that single choice moved an answer by 13.4 kcal/mol and turned an apparent disaster into an agreement. There is no published threshold for CI2. So rather than pick one and present it as the answer, we swept the whole plausible range. Every threshold disagrees. The conclusion does not depend on the choice, which is the only reason we are willing to state it.
And this is not a convergence failure. The run made 167 round trips with a per-replica median of 10, which is the test an earlier run in this project failed with zero. The sampling is sound. The force field is wrong for this molecule.
Why the negative result is the useful one
If CI2 had agreed, this page would be about a method that predicts protein stability from structure. It would have been a more exciting page and a less honest one, because a single success is thin evidence either way.
What we actually have is sharper: the calibration is per protein. You can tune this class of model to reproduce a measurement on one molecule, which is the 1SHF result. Whether the tuned model then predicts any other measurement on that molecule we have not yet tested. What you cannot do is carry the tuning to the next molecule and expect an answer. On the evidence here, each new protein needs its own calibration, and that calibration needs a published value to calibrate against, which is a real limit on using this to predict anything unmeasured.
What it cost
| Quantity | Value |
|---|---|
| Compute spent | 3,817 core-hours |
| Equivalent rented, interruptible | about $70 |
| Machine | one laptop, 18 cores |
| Toolkit directories runnable | 5 of 6, 6 of 6 with substitutes |
| Shipped examples reproduced | 11 |
The headline comparison cost about 160 core-hours at 36,002 exchanges, and about 970 for the whole 216,014-exchange run, where it held at z = −0.06. The expensive part of this work was not the compute. It was finding out which of our own numbers to believe.
What this does not establish
It does not establish that the model is good. It matches one laboratory value in a test it could not have failed, matches one of the toolkit's own values provisionally, and misses a third by a wide margin. The 3CLA melting point comparison has no published uncertainty attached, so no agreement can be tested there at all.
It does not establish that calibrations never transfer. We tested one transfer, between two small two-state folders, and it failed decisively. One failure is not a law, and a different pair might behave differently.
It does not establish why it failed. We can say the model under-stabilises CI2 at this strength. We cannot yet say whether that is the force field, the contact map, the absence of solvent, or something about CI2 in particular.
And the published CI2 value was measured at pH 6.3. Our model has no pH: no titratable groups, no solvent, nothing for a proton to do. We are comparing a pH-dependent measurement against a model that cannot represent pH, and the target temperature of 298.15 K is assumed rather than read from the paper. Neither of those explains a gap this size, but both are real and both are ours.
Where the working is
- The full record, every number traced to the file that produced it.
- Six times we were wrong, the corrections we published on the way here, including the 13.4 kcal/mol one.
- What it takes to run somebody else's research code, the eight defects and two blockers.
- Four walls we hit, what stopped us and what each would take.