[phenixbb] Planarity restraint ignored when dist_esd equal to 0
Pavel Afonine
pafonine at lbl.gov
Tue Sep 27 19:33:54 PDT 2016
Hello Takanori,
thanks for feedback! It is always very much appreciated!
> I have recently encountered a similar situation, where I
> took a CIF file for a ligand from the CCP4 monomer library
> and found planarity restraints are silently ignored in PHENIX.
> This is very confusing and dangerous.
I can see how this may be confusing so I agree with Keitaro that a
warning would be good to have. Also, phenix.refine always outputs .geo
file that lists all the restraints for all atoms, which allows to check
presence or absence of particular restraints for atoms in question.
> > with this logic you interpret the number and not use it. Follow this
> > road you may end up with odd things, like let's treat 0 as
> > "constraints", 0.1 as "strong restraints", 1.5 as "weak restraints",
> etc
> > and in each case do something very different - I hope you get the idea.
>
> The "zero sigma" limit of a Gaussian is an infinitely sharp
> probability distribution (= delta function). So I think the transition
> from restraints to constraints is not "odd" but quite natural both
> mathematically and semantically.
Theoretically it is all right of course, but in terms of code it may be
quite different pieces of code. It is much less cryptic if a decision to
use one option (restraints) or another (constraints) is a made based on
a clearly named parameter rather than inferring from interpreting a number.
> Anyway, if PHENIX does not want to treat 0 as constraints,
> I would prefer an error stop, not just a warning,
It is always a balance between being conservative vs relaxed, automated
vs less automated. There is no perfect solution for drawing the line,
Iafraid. I'm sure if we make phenix.refine stop every time it sees esd=0
there will be a great deal of people flooding mailboxes saying "why it's
stopping, I don't care, I want my refinement job finished", etc.
phenix.refine favors automation paradigm while providing a lot of room
for flexibility. You can define custom bonds, angles, torsions, planes,
groups of planes to be parallel, as many as you wish and for any atoms
of your choice. And you can make sure that what you defined was actually
used by inspecting .geo files.
A practical solution to "esd=0 issue" that is somewhere in between the
two extremes ("stop if esd=0" vs "ignore and keep going") is to replace
esd=0 with say esd=0.01 (or whatever is used for rings in TYR or PHE).
This will let program going as before and you will get your plane that I
bet one wouldn't be able to tell that it is not perfect just by looking
at it!
> because most people
> don't mean "ignore this record" by zero sigma.
Frankly, I don't know and I wouldn't be so sure. I'd say many wouldn't
tell (and care about) the difference between restraints vs constrains,
not to mention nuances about meaning of esd=0, which is fine.
Those who don't mean "ignore this record" are probably trained to use
the program X, where this assumption is made. Perhaps it's not too
unreasonable to realize that not everything adopted in program X is
equally and fully adopted in program Y.
I can make a list of hundreds other similar decisions that are made
under the hood of phenix.refine automatically.
- Do hydrogens participate in bulk-solvent mask calculation?
- Do vdw radii of non-H change when H are used in refinement?
- What program does if it sees non-positive definite anisotropic ADP or
any of T, L or S matrix?
- What happens if input PDB has ANISOU records, TLS is not requested and
resolution is not great for anisotropic refinement?
- Are H treated as riding or freely refined if data resolution is 0.73A?
- What data resolution is used when you request to do rigid-body refinement?
- What refinement target is used (LS or ML) if there is less than 50
test (free) reflections?
... Yep, I can go on for 100+ more items like this. If every time we
stop and ask users what to do that would not go well!
All the best,
Pavel
More information about the phenixbb
mailing list