Inferring class average orientations

The previous section showed how to infer a mixture model given images with known orientations.

This section considers the opposite situation: given a mixture model, how to infer the orientation of each image.

As explained in Chapter 2, theith image orientation is modeled as a rotationR_i. The rotation

1/10

Fourier shell correlation (FSC) ^0∘₁₅∘

30^∘

Fourier shell correlation (FSC) ^0∘

15^∘ 30^∘ 45^∘ 60^∘

0^◦ 15^◦ 30^◦ 45^◦ 60^◦

Figure 3.13: Comparing the effect of the missing cone on the mixture model reconstruction algorithm and on direct Fourier inversion. (Top) Results of the mixture model reconstruction algorithm. (Bottom) Results using direct Fourier inversion. (Right) Each angle denotes the size of the missing cone in Fourier space. A missing cone of 30^◦ corresponds to a maximum tilt angle of 60^◦. (Left) Each reconstruction is compared to a reference GroEL structure at 15 ˚A by computing an FSC curve. Comparing the FSC curves shows that at low resolution (below 20 ˚A) the mixture model algorithm performs better, while at higher resolution direct Fourier inversion is better. At low resolutions, larger missing cones have a stronger effect on the direct Fourier inversion result than the mixture model result. Visually, the mixture model results also appear to be less affected by the missing cone.

describes both the direction in which the mixture model is projected (two parameters), and the in-plane rotation of the projected image (one parameter). Other reconstruction algorithms often treat the projection direction and the in-plane rotation separately.

As in the previous section, the translationsti are fixed: ti= 0 for each i.

Estimating the rotations forms part of many of the reconstruction algorithms described in Chapter 1. For instance, one of the two steps of every projection matching iteration is to update the orientation parameters for each image based on the current estimate of the density map.

This can be either a global or a local update. During the first projection matching iterations, the density map is projected along a grid covering all possible directions, to estimate the globally best rotation for each image. During later iterations, the global grid is replaced by a local one surrounding the current best estimate of the orientation parameters for each image.

0 1 10 50 100 200 300 400 500 1000 input

Figure 3.14: Estimating Pol II image orientation using local Gibbs sampling. The goal is to estimate the orientation of a single image (right-most column). The initial images (left-most column) are projections of the mixture model after applying different random rotations. Even after 1000 Gibbs sampling steps (second-last column), only two of the rotations have converged to the correct solution. The others are stuck in local optima. This shows that the local Gibbs sampler can only find the correct rotation if it is initialised with a nearby rotation.

Similarly, the rotation sampler for mixture models can be either local or global. The local rotation sampler is derived from the Gibbs sampler for the full model (Section 2.4) in the same way that the mixture model sampler (Section 3.2) was derived: by fixing the known parameters and removing their updating steps from the algorithm.

The parameters that remain are the assignments z, the missing components x^m, and the rotationsR. These are sampled in turn from the same conditional distributions as before (Eqns.

2.85, 2.86, 2.92).

Each sampling iteration ofz,x^mandRforms one step of the Gibbs sampler for rotations. The algorithm is initialised by sampling each rotation from the uniform prior over SO(3). Although the Gibbs sampler does not converge to a single rotation, it usually only changes by small amounts once the burn-in period has passed. The algorithm can therefore be evaluated by testing if the last sampled rotation is in the vicinity of the true rotation or not.

Testing the algorithm with Pol II images shows that it does not always converge to the correct rotation (Fig. 3.14). It converges to a local optimum, which coincides with the global optimum only if the initial random rotation happens to be near the true rotation.

One reason why the Gibbs sampler often gets stuck in a local optimum is the following:

After sampling the assignments z, each count x_ijl in the image has been assigned to a specific component, given by z_ijl. The conditional distribution for R_i depends on this assignment, i.e.

it assumes that it is correct. It therefore assigns a very low probability to rotations where the count would be very far away from the projection of its assigned component. But if the current estimate of the rotationR_i is far from the true rotation, then many of the sampled assignments are also wrong, and consequently the true rotation may have a very low likelihood under the conditional distribution for Ri.

This effect of the Gibbs sampler getting stuck in local optima becomes more severe as the number of components increases, and the component size decreases.

To motivate a different approach to rotation sampling, note that as before, the Gibbs sampler is being used to sample from the extended posterior p(θ,Z|D), where the model parameters θ are now just the rotations. As an alternative, consider the posterior that does not involve any latent variables: The posterior forR_i is not a standard distribution from which samples can easily be drawn;

this was the motivation for introducing latent variables in Section 1.4.1. But we can nevertheless sample from Eqn. 3.8 by approximating it with a discrete distribution. The approximation is formed by sampling a large numberNRof rotationsRmuniformly from SO(3), sayNR = 10000.

The values of the posterior of R_i at these rotations form the weights of the discrete approxima-tion:

To sample a rotation Ri from Eqn. 3.9, just sample m from the categorical distribution with weights w_m, and let R_i =R_m be the desired sample.

This approach using the posterior distribution forR_i (Eqn. 3.8) and its discrete approxima-tion (Eqn. 3.9), will be referred to as the global rotaapproxima-tion sampler. It is very similar to the global rotation update step in projection matching.

−1.5 −1.0 −0.5 0.0 0.5 1.0 1.5 Angle error (degrees)

0 2 4 6 8 10 12

Number of rotation estimates ^alphabeta_gamma

0 1 10 50 100 input

A B

Figure 3.15: Combining global and local rotation sampling to estimate the Pol II image orient-ation from Fig. 3.14. (A) The initial rotorient-ations are again sampled randomly (first column). The first step (second column) is a global rotation sampling step, which finds an rotation close to the true rotation. During the remaining local Gibbs sampling steps, each rotation converges to the vicinity of the true rotation. (B) For each of the 25 images, we compare the Euler angles of the true rotation and the estimated rotation. All rotation estimates are very accurate, with most of the angular errors <1^◦.

One difference between the local and global approaches to rotation sampling is that global rotation sampling is much more computationally intensive. Fortunately, one global rotation sampling step is usually enough to find a rotation in the vicinity of the true rotation. The proposed algorithm for estimating rotations is to start with a single global rotation sampling step, followed by multiple local rotation sampling steps to converge to the solution.

Fig. 3.15 shows that this algorithm successfully converges to the correct rotation in all tested examples.

3.4 Inferring mixtures from class averages with unknown