From babff6d8bff8f0f039e510f6b12b09ca381df4ea Mon Sep 17 00:00:00 2001 From: Filip Leonarski Date: Mon, 12 Oct 2026 16:47:07 +0200 Subject: [PATCH] rugnux --af-model: the null keeps the copies' arrangement The null a placement must beat turned the whole assembly about its common centroid, so with several copies it also moved them off the arrangement the search found. Under tNCS, two parallel copies a tNCS vector apart fit the data whatever the model is, and the null threw that fit away: wrong deposited models reached TF Z 14-47 and 8-34 sd over the null on tNCS crystals, and two of them decided the space group (SPACE_GROUP, a reflection file in the wrong group). Each copy is now turned about its own centroid, all by the same rotation, so only the model's own orientation is random. On a control of 38 open-arm sets (deposited chain randomly re-oriented as search model, 3 unrelated deposited models and a fold-scrambled decoy per set) the wrong models' null falls to 1.4-4.2 sd on the tNCS sets (one at 7.2, left UNDECIDED), both wrong decisions become NOT_PLACED, and no right-model or AlphaFold-model decision changes (nulls stay 15-107 sd). The 193-set AF survey replay loses no GREEN (two UNDECIDED become NOT_PLACED); the four twin/tNCS targets stay correct. Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01EizBKhTqrAYago9KJG3dA3 --- .../MolecularReplacement.cpp | 33 ++++++++++++------- 1 file changed, 22 insertions(+), 11 deletions(-) diff --git a/image_analysis/structure_refinement/MolecularReplacement.cpp b/image_analysis/structure_refinement/MolecularReplacement.cpp index 23e1b4600..971f91405 100644 --- a/image_analysis/structure_refinement/MolecularReplacement.cpp +++ b/image_analysis/structure_refinement/MolecularReplacement.cpp @@ -459,11 +459,14 @@ gemmi::Fractional centroid(const gemmi::Structure &st) { // about its centroid, in standard deviations of the latter: the null model validation builds for a // model that claims something (ModelValidation.cpp). An orientation the crystal cannot tell from the // placement's - within the rigid body's reach of one of its symmetry classes - is drawn again. +// With several copies each is turned about its own centroid, all by the same rotation: the copies keep +// their arrangement, and with it the fit that tNCS gives any model placed parallel at its vector, so that +// only the model's own orientation is random. // `translate`: where given, each replicate is placed where the translation function puts the model in its // orientation (the centroids, fractional, of the model in each of the orientations given), as the search // placed the real one, instead of at the real placement's centroid; the search chose that position among // all others for the real orientation, and a null that is not given the same choice reads high. -double null_sigma(const gemmi::Structure &placed, const gemmi::Mat33 &orientation, const gemmi::SpaceGroup &sg, +double null_sigma(const gemmi::Structure &placed, int copies, const gemmi::Mat33 &orientation, const gemmi::SpaceGroup &sg, const std::vector &merged, double d_min, RigidBodyGPUPool *pool, size_t nthreads, ScoreEngines *scorers, const std::function(const std::vector &)> &translate) { @@ -497,15 +500,23 @@ double null_sigma(const gemmi::Structure &placed, const gemmi::Mat33 &orientatio workers.push_back(std::async(std::launch::async, [&] { for (int i = next++; i < NULL_REPLICATES; i = next++) { gemmi::Structure rep = placed; - std::vector pos = ModelPositions(rep.models[0]); - gemmi::Vec3 c; - for (const gemmi::Position &p : pos) - c += p; - c *= 1.0 / static_cast(pos.size()); - const gemmi::Vec3 to = translate ? gemmi::Vec3(placed.cell.orthogonalize(where[i])) : c; - for (gemmi::Position &p : pos) - p = gemmi::Position(draws[i].multiply(gemmi::Vec3(p) - c) + to); - SetModelPositions(rep.models[0], pos); + std::vector &chains = rep.models[0].chains; + for (int k = 0; k < copies; k++) { + // Copy k: the chains SplitCopies gives it. + std::vector atoms; + for (size_t ch = 0; ch < chains.size(); ch++) + if (static_cast(ch * copies / chains.size()) == k) + for (gemmi::Residue &r : chains[ch].residues) + for (gemmi::Atom &a : r.atoms) + atoms.push_back(&a); + gemmi::Vec3 c; + for (const gemmi::Atom *a : atoms) + c += a->pos; + c *= 1.0 / static_cast(atoms.size()); + const gemmi::Vec3 to = translate ? gemmi::Vec3(placed.cell.orthogonalize(where[i])) : c; + for (gemmi::Atom *a : atoms) + a->pos = gemmi::Position(draws[i].multiply(gemmi::Vec3(a->pos) - c) + to); + } r_work[i] = refine_placement(rep, sg, coarse, null_d_min, pool, nthreads, scorers).r_work; } })); @@ -2440,7 +2451,7 @@ MolecularReplacementResult search_model(const SearchModel &model, const std::vec out.push_back(p.t); return out; }; - c.null_sigma = null_sigma(*c.placed, orientations[best], c.space_group, *merges[merge].merged, + c.null_sigma = null_sigma(*c.placed, c.copies, orientations[best], c.space_group, *merges[merge].merged, result.d_min_score, engines.get(), nthreads, &scorers, translate); result.t_null += seconds_since(t_null); if (c.null_sigma >= NULL_SIGMA_MIN)