In helical reconstruction, psi prior angles are not flipped 180 degrees when needed

[thurberk]

In short, in helical analysis, when a helical structure has a direction opposite to the initial model, Relion is not flipping the AnglePsiPrior by 180 degrees, in a Refine3D job.

This has bad consequences if local angular averaging is triggered either by the initial parameters, or during the Relion Refine3D job.

In more detail:
For helical analysis, when Relion is doing non-local angular averaging, it checks Psi angles within the range from AnglePsiPrior, and also within the range from AnglePsiPrior + 180 degrees.
When Relion is doing local angular averaging, it only checks Psi angles within the range from AnglePsiPrior. So, if the AnglePsiPrior is not flipped by 180 degrees for helices that are going the opposite direction, when local averaging is triggered only the wrong direction of the Psi angle is checked.
Currently (in Relion 5 and probably earlier versions) the AnglePsiPrior is not being flipped by 180 degrees, even when most particles in the helix have preferred orientation of Psi ~180 degrees away from the prior.

There is also a closely related issue where the PsiFlipRatio variable is not being updated if the Psi angle is not flipped. In this case, the PsiFlipRatio variable stays at 0.5, for helices that are going the same direction as the model.

Before Relion is updated, the main bad effects can be avoided by avoiding local angular searches. To be explicit, this means: On the Auto-sampling tab of a 3D run, the “Local searches from auto-sampling” should be set small enough that it will not be reached during the run.

I think the fix in the code is thankfully pretty simple as detailed below:

helix.cpp: Lines 4136 to 4145
// Change the polarity of the entire helix if psi_flip_ratio is larger than 0.5
if (psi_flip_ratio > 0.5)
{
for (int id = sid; id <= eid; id++)
{
flipPsiTiltForHelicalSegment(list[id].psi_prior_deg, list[id].tilt_prior_deg,
list[id].psi_prior_deg, list[id].tilt_prior_deg);
list[id].psi_flip_ratio = (1. - psi_flip_ratio);
}
}

Should be:
// Change the polarity of the entire helix if psi_flip_ratio is larger than 0.5
if (psi_flip_ratio > 0.5)
{
for (int id = sid; id <= eid; id++)
{
flipPsiTiltForHelicalSegment(list[id].psi_prior_deg, list[id].tilt_prior_deg,
list[id].psi_prior_deg, list[id].tilt_prior_deg);
flipPsiTiltForHelicalSegment(list[id].psi_deg, list[id].tilt_deg,
list[id].psi_deg, list[id].tilt_deg);
list[id].psi_flip_ratio = (1. - psi_flip_ratio);
}
}
else
{
for (int id = sid; id <= eid; id++)
{
list[id].psi_flip_ratio = psi_flip_ratio;
}
}

[biochem-fan]

Thank you very much. Because my experience with helical refinement is limited, I assigned @scheres to this issue.

I feel (but not certain) that your second point about PsiFlipRatio is right. On the other hand, isn't the first point more of a semantic issue? I guess local search assumes that angles and priors in the input STAR file are more-or-less correct and that there is no need to search the other side (+ 180).

Just clarification: do you mean local angular search when you mention "local angular averaging"?

[thurberk]

biochem-fan
Thanks for your note.
Yes, you are correct. I should have written "local angular search," not "local angular averaging."

I don't think it is just a semantic issue.
One way the problem can occur is if Relion switches to local angular search in the middle of the Refine3D run. That happens when Relion hits the value in the input "Local searches from auto-sampling" on the Auto-sampling tab in a Refine3D run. At that point, Relion stops checking Psi angles near PsiPrior + 180, and this leads to incorrect values for helices/fibrils that were picked "backwards" from the model.
I think users could run into this problem without realizing it, and they don't have control over the star file in the middle of the run.

Alternatively, users can ask for local angular search from the beginning of the run by setting the "Local searches from auto-sampling" value above the "Initial angular sampling" value. In this case, the user could edit the star file PsiPriors from a previous run outside of Relion, but I still don't think it is a good idea to expect users to realize they need to.

The counter-argument that I see is that the PsiPrior values should always stay at the initial value generated by the picking direction, but I see that as less important. Knowing whether the Psi angle has been flipped by ~180 degrees from the initial picking does have value, because it tells you what the sign of the Rot angle dependence should be moving down the helix (to state that another way, whether the Rot angle should be increasing or decreasing from one particle to the next in a slowly twisting helix.) In one of my own modified versions of Relion, I created a PsiPriorStart variable to hold the initial value of PsiPrior, then Relion can tell whether the PsiPrior has been flipped.

[thurberk]

So, I would advocate for a PsiPriorStart variable, which can hold the initial PsiPrior from picking.
That way you can get both advantages:
The PsiPrior can be flipped 180 degrees if needed without losing the information of the initial helix direction.
Then, local angular searches can be done with the Psi angle.
And the PsiPriorStart can be compared to the current PsiPrior or Psi angle to determine what the relative sign of Rot angle dependence should be down the helix/fibril.

[biochem-fan]

One way the problem can occur is if Relion switches to local angular search in the middle of the Refine3D run. That happens when Relion hits the value in the input "Local searches from auto-sampling" on the Auto-sampling tab in a Refine3D run. At that point, Relion stops checking Psi angles near PsiPrior + 180, and this leads to incorrect values for helices/fibrils that were picked "backwards" from the model.

Isn't PsiPrior already in the correct direction when a refinement job switches from global to local search?

The code you cited seems to update (flip) PsiPrior. Isn't this triggered during global search?
(I am not very familiar with helical refinement code in RELION)

relion/src/helix.cpp

Lines 4136 to 4145 in adfec82

	// Change the polarity of the entire helix if psi_flip_ratio is larger than 0.5
	if (psi_flip_ratio > 0.5)
	{
	for (int id = sid; id <= eid; id++)
	{
	flipPsiTiltForHelicalSegment(list[id].psi_prior_deg, list[id].tilt_prior_deg,
	list[id].psi_prior_deg, list[id].tilt_prior_deg);
	list[id].psi_flip_ratio = (1. - psi_flip_ratio);
	}
	}

[thurberk]

Isn't PsiPrior already in the correct direction when a refinement job switches from global to local search?

Unfortunately, PsiPrior is not flipped to the correct direction for "backwards" helices, and this causes the main problem.

The code you cited does flip the PsiPrior, but this comes at the beginning of the updatePriorsForOneHelicalTube function. If you are doing angular averaging of the Psi angle over the helix, the averaging of the Psi angle occurs later in this function and uses the Psi angle which is not flipped. So, the PsiPrior value is overwritten later in the function by an averaged unflipped Psi value.

Another relevant point is that all the Psi values at this point in the function are unflipped values, because any flipped Psi values are flipped back by this piece of code which comes in this function before both the section you quoted and the angular averaging.

relion/src/helix.cpp

Lines 4117 to 4121 in adfec82

	if (list[id].has_wrong_polarity)
	{
	flipPsiTiltForHelicalSegment(list[id].psi_deg, list[id].tilt_deg, list[id].psi_deg, list[id].tilt_deg);
	nr_opposite_polarity++;
	}

So, in summary, the code section you cite should flip both the PsiPrior and the Psi angle, so if the Psi angles are averaged later in the function and then put in PsiPrior, that new PsiPrior value is also a flipped value.

[biochem-fan]

Okay. I understood your point.

@scheres, this arguments sound valid to me. Do you agree merging the proposed fix?

[scheres]

Sorry for the delay in dealing with this. The past 2 months have been busy. Kent, thank you very much for going through the code so carefully. Your careful contributions are much appreciated.

Am I correct that this bug does NOT affect normal local angular searches (in single-particle or helical refinement)? This only becomes a problem when local averaging of orientations is switched on on the Helix tab, right? The discussion above seems to suggest otherwise, but then I don't understand the problem.

Anyway, I am testing your suggested fix. More soon!

[thurberk]

Am I correct that this bug does NOT affect normal local angular searches (in single-particle or helical refinement)? This only becomes a problem when local averaging of orientations is switched on on the Helix tab, right?

Yes, I agree that this should only be a problem when local averaging is switched on on the Helix tab.
The problem is: When the local averaging is on, the (unflipped) Psi angles are used in the averaging to generate a new (unflipped) PsiPrior value, which replaces the PsiPrior value from earlier in the function.
If the local averaging doesn't occur, the earlier (flipped) PsiPrior value will stay.