Post-refinement (xsphere) feedback and questions

[biochem-fan]

This is not a bug report but feedback and questions.

I am getting a very promising result using partialator (xsphere model). The resolution as judged by CC1/2 ~ 0.5 improved from 2.45 to 2.20 Å for a real dataset (not lysozyme).

The dataset was integrated with the default bandwidth (1E-8). So I set --force-bandwidth 0.0046 to partialator to provide a more realistic initial value. I noticed that the convergence was slower without this (5 vs 10 cycles). Nonetheless, after many cycles (10 cycles), both jobs with/without --force-bandwidth 0.0046 reached similar CC1/2, Rsplit and I/sig. This probably means that the algorithm can eventually find the right bandwidth even when starting from 1E-8. Interestingly, however, the multiplicity without --force-bandwidth 0.0046 was only about half of that with --force-bandwidth 0.0046. Apparently more reflections were rejected with the narrower initia bandwiwdth (perhaps due to too small partialities?). Puzzlingly, the number of rejected crystals were larger with the larger initial bandwidth (2214 B too big for 0.0046 vs 843 B too big for 1E-8). I don't know how to interpret this. Perhaps more reflections that originated from the tails of the spectrum could be incorporated with a broader initial bandwidth, but they were not properly modeled and not contributing useful signal to the final dataset? Do you have any thoughts?

Of course, these metrics merely reflect internal consistency, so I have to check atomic model refinement as well.

More generally, I would like to understand which parameter affects what. The bandwidth seems to be used in PinkIndexer and partialator, but not by profile radius calculation before integration. Does it affect prediction? I see the new spectrum object is used in the code but am not sure if the work is complete. The divergence is used by PinkIndexer, but not by partialator (although it is required), profile radius calculation before integration or prediction. The mosaicity does not seem to be used at all (!). Am I correct?

Another topic: I remember that partialator used to refine not only crystal directions but also unit cell parameters in earlier versions, but this feature was dropped. Is this because it turned out to be harmful (i.e. destabilized refinement)? Now that we have a mechanism to respect Bravais lattice constraints, it might be interesting to revisit this.

[taw10]

Post-refinement in 0.12.0 (back to at least 0.10.0) is only refining four parameters: crystal tilt (two directions), profile radius and wavelength. Bandwidth doesn't get refined at all, so it's interesting that it refines to similar self-consistency.

I'm not sure why there would be lower multiplicity with a larger bandwidth. I don't think it's because of the partiality model itself, because that should go in the opposite direction. Can the number of rejected crystals (not reflections) explain the reduction in multiplicity? I think that "B too big" is often a symptom rather than a cause. Something was "wrong" with those crystals, and the scaling calculation converged to a silly value.

These days I'm not so sure that the convergence properties of post-refinement are that good. The current way used by all software, alternately refining and merging a new reference, seems like it should work. But similar procedures in at least two other areas (MIR heavy atom refinement and detector geometry refinement) turned out to be convergence disaster areas even though it looked like they were working. More discussion about this in my upcoming paper.

As for which parameter affects what, it depends on the partiality model. Spot prediction uses the 'ggpm' partiality model. Prediction parameter determination only fits the profile radius, not mosaicity nor divergence nor bandwidth.

ggpm uses bandwidth via spectrum object (including multiple Gaussians, if provided), profile radius, but not mosaicity nor divergence.
xsphere uses bandwidth via spectrum object (can deal with an arbitrary spectrum), profile radius, mosaicity, but not divergence.
offset uses profile radius and mosaicity, but not spectrum (it's a completely monochromatic model).
unity ignores everything, sets all partialities to 1.

Mosaicity could be added to ggpm quite easily (two lines, in fact). And looking through the code, I don't think divergence is used anywhere at all.

I definitely agree about refining cell parameters again. I think it would be much more stable now, and a lot of the refinement machinery can be re-used.

[biochem-fan]

Thank you very much for your feedback.

Bandwidth doesn't get refined at all, so it's interesting that it refines to similar self-consistency.

Oops. I somehow thought the bandwidth was refined. Probably I confused with early days of Helen's cppxfel.

I'm not sure why there would be lower multiplicity with a larger bandwidth.

No, "the multiplicity without --force-bandwidth 0.0046 was only about half of that with --force-bandwidth 0.0046", so the multiplicity was larger with a larger bandwidth. Your intuition is correct.

These days I'm not so sure that the convergence properties of post-refinement are that good. The current way used by all software, alternately refining and merging a new reference, seems like it should work. But similar procedures in at least two other areas (MIR heavy atom refinement and detector geometry refinement) turned out to be convergence disaster areas even though it looked like they were working. More discussion about this in my upcoming paper.

This is an interesting observation and I am looking forward to reading your paper.

It might be useful to introduce stochasticity to the optimizer. The initial model problem in CryoEM SPA has exactly the same mathematical structure as post-refinement and various stochastic optimization algorithms solved this long-standing problems. We can try subset EM (do not update paramaters of all crystals but only a subset in an iteration) and/or stochastic hill climbing (do not always choose the best direction). Although computationally more demanding, we can also try full marginalization over crystal parameters or do stochastic EM (sample crystal parameters from its posterior distributions).

As for which parameter affects what, it depends on the partiality model. Spot prediction uses the 'ggpm' partiality model. Prediction parameter determination only fits the profile radius, not mosaicity nor divergence nor bandwidth.

Thank you very much for clarification.

I definitely agree about refining cell parameters again. I think it would be much more stable now, and a lot of the refinement machinery can be re-used.

Do you have time to work on it? After my beam time next week, I myself might be able to hack the code to refine a, b, c and see if it is worthwhile (my crystals are orthorhombic).

[taw10]

I'm planning to work (already am) on post-refinement, or more widely crystal modelling, but I had a different approach in mind than adding parameters to the current algorithm. So, you should feel free to try adding them and see what happens.

[biochem-fan]

A short update:

I tested refinement of unit cell parameters (no angles supported at the moment). Apparently they are strongly correlated with crystal orientations. If I refined the orientation first, the lengths did not change. If I refined the lengths first, the orientation hardly moved. Probably I need a better optimiser than simple step searches. Unfortunately I did not have time to go deeper and suspended my effort. This turned out non-trivial.