Common lines ab-initio reconstruction of $D_2$-symmetric molecules

TL;DR

This paper introduces a novel common lines algorithm specifically designed for ab-initio reconstruction of molecules with D2 symmetry in cryo-electron microscopy, leveraging all pairwise and within-image common lines.

Contribution

The work presents the first common lines algorithm tailored for D2-symmetric molecules, extending existing methods from cyclic to dihedral symmetry.

Findings

Successfully applied to experimental cryo-EM data

Demonstrates accurate structure determination for D2-symmetric molecules

Extends common lines methodology to new symmetry class

Abstract

Cryo-electron microscopy is a state-of-the-art method for determining high-resolution three-dimensional models of molecules, from their two-dimensional projection images taken by an electron microscope. A crucial step in this method is to determine a low-resolution model of the molecule using only the given projection images, without using any three-dimensional information, such as an assumed reference model. For molecules without symmetry, this is often done by exploiting common lines between pairs of images. Common lines algorithms have been recently devised for molecules with cyclic symmetry, but no such algorithms exist for molecules with dihedral symmetry. In this work, we present a common lines algorithm for determining the structure of molecules with $D_{2}$ symmetry. The algorithm exploits the common lines between all pairs of images simultaneously, as well as common lines within each image. We demonstrate the applicability of our algorithm using experimental cryo-electron microscopy data.

Equations

$$R_{i}=\left(\begin{array}[]{ccc}|&|&|\\ R_{i}^{(1)}&R_{i}^{(2)}&R_{i}^{(3)}\\ |&|&|\end{array}\right)\in SO(3),\quad i\in[N],$$

$$P_{R_i}(x,y)={\int }_{-\infty }^{\infty }\varphi (R_{i}r)dz={\int }_{-\infty }^{\infty }\varphi (x{R}_i^{(1)}+y{R}_i^{(2)}+z{R}_i^{(3)})dz.$$

$$g_2=\left(\begin{array}{ccc}[r]1& 0& 0\\ 0& -1& 0\\ 0& 0& -1\end{array}\right),g_3=\left(\begin{array}{ccc}[r]-1& 0& 0\\ 0& 1& 0\\ 0& 0& -1\end{array}\right),g_4=\left(\begin{array}{ccc}[r]-1& 0& 0\\ 0& -1& 0\\ 0& 0& 1\end{array}\right)$$

$$\varphi (r)=\varphi (g_{1}r)=\varphi (g_{2}r)=\varphi (g_{3}r)=\varphi (g_{4}r)$$

$$P_{R_i}(x,y)={\int }_{-\infty }^{\infty }\varphi (R_{i}r)dz={\int }_{-\infty }^{\infty }\varphi (g_m{R}_{i}r)dz=P_{g_m{R}_i}(x,y)$$

$$P_{R_i}(x,y)={\int }_{-\infty }^{\infty }\varphi (R_{i}r)dz={\int }_{-\infty }^{\infty }\psi (J{R}_{i}r)dz={\int }_{-\infty }^{\infty }\psi (J{R}_{i}JJr)dz,$$

$${\int }_{-\infty }^{\infty }\psi ((J{R}_{i}J)Jr)dz={\int }_{-\infty }^{\infty }\psi (J{R}_{i}J\left(\begin{array}{c}x\\ y\\ z^{\prime }\end{array}\right))d{z}^{\prime }={\tilde{P}}_{J{R}_{i}J}(x,y),$$

$${\hat{P}}_{R_i}({\omega }_x,{\omega }_y)=\hat{\varphi }({\omega }_x{R}_i^{(1)}+{\omega }_y{R}_i^{(2)}),({\omega }_x,{\omega }_y)\in {\mathbb{R}}^2.$$

$$q_{ij}=\frac{R_i^{(3)}\times R_j^{(3)}}{\Vert R_i^{(3)}\times R_j^{(3)}\Vert }.$$

$$q_{ij}=\cos {\alpha }_{ij}{R}_i^{(1)}+\sin {\alpha }_{ij}{R}_i^{(2)}=\cos {\alpha }_{ji}{R}_j^{(1)}+\sin {\alpha }_{ji}{R}_j^{(2)},$$

$${\hat{P}}_{R_i}(\xi \cos {\alpha }_{ij},\xi \sin {\alpha }_{ij})={\hat{P}}_{R_j}(\xi \cos {\alpha }_{ji},\xi \sin {\alpha }_{ji}),\xi \in \mathbb{R}.$$

$$q_{ij}^m=\frac{R_i^{(3)}\times g_m{R}_j^{(3)}}{\Vert R_i^{(3)}\times g_m{R}_j^{(3)}\Vert },m=1,2,3,4,$$

$$\begin{array}{rl}C(R_i,g_m{R}_j)& =\left(\begin{array}{c}\cos {\alpha }_{ij}^m,\sin {\alpha }_{ij}^m\end{array}\right)=\left(\begin{array}{c}⟨R_i^{(1)},q_{ij}^m⟩,⟨R_i^{(2)},q_{ij}^m⟩\end{array}\right),\\ C(g_m{R}_j,R_i)& =\left(\begin{array}{c}\cos {\alpha }_{ji}^m,\sin {\alpha }_{ji}^m\end{array}\right)=\left(\begin{array}{c}⟨g_m{R}_j^{(1)},q_{ij}^m⟩,⟨g_m{R}_j^{(2)},q_{ij}^m⟩\end{array}\right),\end{array}$$

$${\hat{P}}_{R_i}(\xi \cos {\alpha }_{ij}^m,\xi \sin {\alpha }_{ij}^m)={\hat{P}}_{R_j}(\xi \cos {\alpha }_{ji}^m,\xi \sin {\alpha }_{ji}^m),\xi \in \mathbb{R},$$

$$q_{ii}^m=\frac{R_i^{(3)}\times g_m{R}_i^{(3)}}{\Vert R_i^{(3)}\times g_m{R}_i^{(3)}\Vert },m=2,3,4.$$

$${\hat{P}}_{R_i}(\xi \cos {\alpha }_{ii}^{1m},\xi \sin {\alpha }_{ii}^{1m})={\hat{P}}_{R_i}(\xi \cos {\alpha }_{ii}^{m1},\xi \sin {\alpha }_{ii}^{m1}),\xi \in \mathbb{R},$$

$$\begin{array}{rl}C(R_i,g_m{R}_i)& =\left(\begin{array}{c}\cos {\alpha }_{ii}^{1m},\sin {\alpha }_{ii}^{1m}\end{array}\right)=\left(\begin{array}{c}⟨R_i^{(1)},q_{ii}^m⟩,⟨R_i^{(2)},q_{ii}^m⟩\end{array}\right),\\ C(g_m{R}_i,R_i)& =\left(\begin{array}{c}\cos {\alpha }_{ii}^{m1},\sin {\alpha }_{ii}^{m1}\end{array}\right)=\left(\begin{array}{c}⟨g_m{R}_i^{(1)},q_{ii}^m⟩,⟨g_m{R}_i^{(2)},q_{ii}^m⟩\end{array}\right),\end{array}$$

$$M_{ij}=R_i^T{R}_j,i,j\in [N].$$

$$M=U^{T}U,$$

$${\mathcal{D}}_c=\{\{Q_l^T{g}_m{Q}_r{\}}_{m=1}^4|Q_l,Q_r\in SO(3),|<Q_l^{(3)},Q_r^{(3)}>|\ne 1\},$$

$${\tilde{q}}_{lr}^m=\frac{Q_l^{(3)}\times g_m{Q}_r^{(3)}}{\Vert Q_l^{(3)}\times g_m{Q}_r^{(3)}\Vert },m\in \{1,2,3,4\},$$

$$C(Q_l,g_m{Q}_r)=\left(\begin{array}{c}\cos {\tilde{\alpha }}_{lr}^m,\sin {\tilde{\alpha }}_{lr}^m\end{array}\right),C(g_m{Q}_r,Q_l)=\left(\begin{array}{c}\cos {\tilde{\alpha }}_{rl}^m,\sin {\tilde{\alpha }}_{rl}^m\end{array}\right),$$

$${\nu }_{n,\theta }(\xi )={\hat{P}}_{R_n}(\xi \cos \theta ,\xi \sin \theta ),\xi \in (0,\infty ),$$

$${\rho }_{ij}({\tilde{\alpha }}_{lr}^m,{\tilde{\alpha }}_{rl}^m)=\frac{{\int }_0^{\infty }({\nu }_{i,{\tilde{\alpha }}_{lr}^m}(\xi ){)}^{\ast }{\nu }_{j,{\tilde{\alpha }}_{rl}^m}(\xi )d\xi }{\Vert {\nu }_{i,{\tilde{\alpha }}_{lr}^m}(\xi ){\Vert }_{L_2}\Vert {\nu }_{j,{\tilde{\alpha }}_{rl}^m}(\xi ){\Vert }_{L_2}},m=1,2,3,4,$$

$${\pi }_{ij}(Q_l,Q_r)=\prod _{m=1}^4{\rho }_{ij}({\tilde{\alpha }}_{lr}^m,{\tilde{\alpha }}_{rl}^m).$$

$${\tilde{q}}_{ll}^m=\frac{Q_l^{(3)}\times g_m{Q}_l^{(3)}}{\Vert Q_l^{(3)}\times g_m{Q}_l^{(3)}\Vert },{\tilde{q}}_{rr}^m=\frac{Q_r^{(3)}\times g_m{Q}_r^{(3)}}{\Vert Q_r^{(3)}\times g_m{Q}_r^{(3)}\Vert },m=2,3,4,$$

$$C(Q_l,g_m{Q}_l)$$

$$C(Q_r,g_m{Q}_r)$$

$${\rho }_{ii}({\tilde{\alpha }}_{ll}^{1m},{\tilde{\alpha }}_{ll}^{m1})$$

$${\rho }_{jj}({\tilde{\alpha }}_{rr}^{1m},{\tilde{\alpha }}_{rr}^{m1})$$

$${\pi }_{ij}(Q_l,Q_r)=\prod _{m=1}^4{\rho }_{ij}({\tilde{\alpha }}_{lr}^m,{\tilde{\alpha }}_{rl}^m)\prod _{m=2}^4{\rho }_{ii}({\tilde{\alpha }}_{ll}^{1m},{\tilde{\alpha }}_{ll}^{m1}){\rho }_{jj}({\tilde{\alpha }}_{rr}^{1m},{\tilde{\alpha }}_{rr}^{m1}).$$

Full text

In this work, we present a common lines algorithm for determining the structure of molecules with $D_{2}$ symmetry. The algorithm exploits the common lines between all pairs of images simultaneously, as well as common lines within each image. We demonstrate the applicability of our algorithm using experimental cryo-electron microscopy data.

1 Introduction

Cryo-electron microscopy (cryo-EM) is a technique for acquiring two-dimensional projection images of biological macromolecules [5]. In this technique, a large number of copies of the same molecule is rapidly frozen in a thin layer of ice, fixing each molecule in some random unknown orientation. The frozen specimen is then imaged by an electron microscope, producing a set of two-dimensional projection images (defined below). Once the imaging orientations of the frozen molecules which produced the images are obtained, the three-dimensional structure of the molecule can be recovered from the projection images by standard tomographic procedures.

Formally, choosing some arbitrary fixed coordinate system of $\mathbb{R}^{3}$, the orientations of the imaged molecules at the moment of freezing can be described by a set of rotation matrices

[TABLE]

where we denote by $[N]$ the set $\{1,\ldots,N\}$, and $SO(3)$ is the group of $3\times 3$ rotation matrices. We denote the density function of the molecule by $\phi(r):\mathbb{R}^{3}\to\mathbb{R}$, where $r=(x,y,z)^{T}$, and by $P_{R_{i}}$ the image generated by the microscope when imaging a copy of $\phi$ rotated by $R_{i}$. The image $P_{R_{i}}$ is then given by the line integrals of $\phi(r)$ along the lines parallel to $R_{i}^{(3)}$, namely

[TABLE]

The orthogonal unit vectors $R^{(1)}_{i}$ and $R^{(2)}_{i}$, which span the plane perpendicular to $R_{i}^{(3)}$, form the $(x,y)$ coordinate system for the image $P_{R_{i}}$ from the point of view of an observer looking from the direction of the electron beam. We refer to $R_{i}^{(3)}$ and its perpendicular plane as the beaming direction and the projection plane of $P_{R_{i}}$, respectively.

We can now state the “orientation assignment problem” as the task of finding a set of $N$ matrices $\{R_{1},\ldots,R_{N}\}\in SO(3)$ such that (1.2) is satisfied for all $i\in[N]$, given only the images $P_{R_{1}},\ldots,P_{R_{N}}$ (in particular, $\phi$ in (1.2) is unknown). In this work, we address the task of determining the orientations of a set of projection images obtained from a molecule with $D_{2}$ symmetry.

In plain language, a $D_{2}$-symmetric molecule has three mutually perpendicular symmetry axes, where after we rotate the molecule by $180^{\circ}$ about any of these axes, the molecule looks exactly the same. Formally, let

[TABLE]

denote the three rotation matrices by $180^{\circ}$ about the $x,y$ and $z$ axes, respectively. For notational convenience we also denote the $3\times 3$ identity matrix by $g_{1}$. Choosing a coordinate system in which the rotational symmetry axes of the molecule coincide with the $x,y$ and $z$ axes, the $D_{2}$ symmetry property implies that

[TABLE]

for any $r\in\mathbb{R}^{3}$. Considering any projection image $P_{R_{i}}(x,y)$, by (1.4) we have

[TABLE]

for $m=2,3,4$. Equation (1.5) shows that a $D_{2}$-symmetric molecule induces an ambiguity in which all orientation assignments of the form $\{g_{m_{i}}R_{i}\}_{i=1}^{N}$, $g_{m_{i}}\in\{g_{1},g_{2},g_{3},g_{4}\}$ are consistent with the same set of images $\{P_{R_{i}}\}_{i=1}^{N}$.

An additional ambiguity inherent to cryo-EM arises from the well known fact that the handedness (chirality) of a molecule cannot be established from its projection images. Denoting by $J=\text{diag}(1,1,-1)$ the reflection matrix through the $xy$-plane, we define by $\psi(r)=\phi(Jr)$ the mirror image of the molecule $\phi(r)$. Since $J^{2}=I$, we also have that $\phi(r)=\psi(Jr)$, and thus, by (1.2) we have

[TABLE]

where $r=(x,y,z)^{T}$. Noting that $Jr=(x,y,-z)^{T}$, and changing the variable $z$ to $z^{\prime}=-z$ we have

[TABLE]

which shows that any projection image $P_{R_{i}}$ of the molecule $\phi(r)$ is identical to the projection $\widetilde{P}_{JR_{i}J}$ of its mirror image molecule $\psi(r)$. Thus, both orientation assignments $\{R_{k}\}_{k=1}^{N}$ and $\{JR_{k}J\}_{k=1}^{N}$ are consistent with the same set of projection images $\{P_{R_{1}},\ldots,P_{R_{N}}\}$.

The “orientation assignment problem” for a $D_{2}$-symmetric molecule can now be stated as the task of finding one of the sets of matrices $\{{R}_{i}\}_{i=1}^{N}$ or $\{J{R}_{i}J\}_{i=1}^{N}$ satisfying (1.2), where each matrix $R_{i}$ can be independently replaced by a matrix $\widetilde{R}_{i}\in\{g_{m}R_{i}\}_{m=2}^{4}$.

2 Common lines and their $D_{2}$ induced geometry

One of the principal approaches to solving the orientation assignment problem, which we also employ in this work, relies on the well know Fourier slice theorem [7], which states that the two-dimensional Fourier transform of a projection image $P_{R_{i}}$ is a central planar slice of the three-dimensional Fourier transform of the molecule $\phi(r)$. Formally, denoting by $\hat{\phi}$ the three-dimensional Fourier transform of $\phi(r)$, and by $\hat{P}_{R_{i}}$ the two-dimensional Fourier transform of the image $P_{R_{i}}$, we have that

[TABLE]

Note that the central slice of $\hat{\phi}$ in (2.1) is spanned by the vectors $R_{i}^{(1)}$ and $R_{i}^{(2)}$ (the first two columns of the rotation matrix $R_{i}$ of (1.1)). Since any two non-coinciding central planes intersect at a single line through the origin, it follows that any pair of transformed projection images $\hat{P}_{R_{i}}$ and $\hat{P}_{R_{j}}$ for which $|\langle R_{i}^{(3)},R_{j}^{(3)}\rangle|\neq 1$ share a unique pair of identical central lines. Formally, consider the unit vector

[TABLE]

By the definition of the cross product, $q_{ij}$ is perpendicular to the normal vectors of both projection planes of $\hat{P}_{R_{i}}$ and $\hat{P}_{R_{j}}$, and thus it gives the direction of their common line (see [11] for a detailed discussion). We can express $q_{ij}$ using its local coordinates on both projection planes by

[TABLE]

where $\alpha_{ij}$ and $\alpha_{ji}$ are the angles between $q_{ij}$ and the local $x$-axes of the planes. Using this notation, the common line property implies that

[TABLE]

Now, let $P_{R_{i}}$ and $P_{R_{j}}$ be a pair of images of a $D_{2}$-symmetric molecule. By (1.5), the image $P_{R_{j}}$ (also denoted $P_{g_{1}R_{j}}$) is identical to the three images $P_{g_{2}R_{j}},P_{g_{3}R_{j}}$ and $P_{g_{4}R_{j}}$. However, each of these four images corresponds to a plane, and by (1.2) the four planes are different from each other since their rotation matrices are different. Each of these planes has a unique common line with the projection plane of ${P}_{R_{i}}$, and thus, we conclude that there are four different pairs of identical central lines between the images $\hat{P}_{R_{i}}$ and $\hat{P}_{R_{j}}$. The directions of all four common lines are given by the unit vectors

[TABLE]

where for convenience we denote $q_{ij}^{1}=q_{ij}$, since $g_{1}$ is the identity matrix. We write the local coordinates of the common lines of the pairs of images $\hat{P}_{R_{i}}$ and $\hat{P}_{g_{m}R_{j}}$, on the respective projection plane of each image (see (2.3)), as

[TABLE]

for $m=1,2,3,4$. Thus, the common line property for $D_{2}$-symmetric molecules implies that

[TABLE]

for $m=1,2,3,4$. Throughout the following sections, we refer to all four common lines as the common lines of the images $P_{R_{i}}$ and $P_{R_{j}}$.

An additional feature of symmetric molecules, and in particular of $D_{2}$-symmetric molecules, is self common lines. For any $i\in[N]$, the images $P_{R_{i}},P_{g_{2}R_{i}},P_{g_{3}R_{i}}$ and $P_{g_{4}R_{i}}$ of a $D_{2}$-symmetric molecule, are identical. However, each image corresponds to a different projection plane, and thus, the projection plane of $P_{R_{i}}$ has a unique common line with each of the planes of $P_{g_{2}R_{i}},P_{g_{3}R_{i}}$ and $P_{g_{4}R_{i}}$. The directions of these common lines are given by

[TABLE]

Thus, there are three different pairs of central lines in $\hat{P}_{R_{i}}$, corresponding to the common lines between the pairs $\{P_{R_{i}},P_{g_{2}R_{i}}\},\{P_{R_{i}},P_{g_{3}R_{i}}\}$ and $\{P_{R_{i}},P_{g_{4}R_{i}}\}$, which satisfy

[TABLE]

for $m=2,3,4$, where

[TABLE]

for $m=2,3,4$. We refer to all three common lines in (2.9) as the self common lines of the image $P_{R_{i}}$.

3 Related work

Many of the current common line methods for orientation assignment are based on the method of angular reconstitution by Van Heel [13]. The core idea of the angular reconstitution method is rooted at the observation that the intersection of any three non-coinciding central planes establishes the orientations of all three, relative to each other. In particular, Van Heel shows how given a triplet of projection images $\{P_{R_{i}},P_{R_{j}},R_{R_{k}}\}$, one can obtain either of the sets of relative rotation matrices $\{R_{i}^{T}R_{j},R_{i}^{T}R_{k},R_{j}^{T}R_{k}\}$ or $\{JR_{i}^{T}R_{j}J,JR_{i}^{T}R_{k}J,JR_{j}^{T}R_{k}J\}$, by using the common lines between $P_{R_{i}},P_{R_{j}}$ and $P_{R_{k}}$. Both choices of relative rotation matrices are equally consistent with the images $P_{R_{i}},P_{R_{j}}$ and $P_{R_{k}}$, and are just the manifestation of the handedness ambiguity discussed in the previous section. The angular reconstitution method then makes the assumption that (without loss of generality) $R_{i}=I$, which immediately establishes $R_{j}$ and $R_{k}$ from $R_{i}^{T}R_{j}$ and $R_{i}^{T}R_{k}$. The orientations of the rest of the images $P_{R_{l}}$ for $l\neq i,j,k$ are obtained by fixing the pair of images $P_{R_{i}}$ and $P_{R_{j}}$, and applying the same method sequentially to each triplet of images $\{P_{R_{i}},P_{R_{j}},P_{R_{l}}\}$ to retrieve $R_{i}^{T}R_{l}$, which immediately establishes $R_{l}$ by $R_{l}=R_{i}^{T}R_{l}$. Note that since all relative rotations and subsequently the rotations themselves were obtained by combining each of the images $R_{l}$ with the same pair of images $\{P_{R_{i}},P_{R_{j}}\}$, the angular reconstitution method ensures that we obtain a hand-consistent assignment, i.e. we either obtain $\{R_{i}\}_{i=1}^{N}$ or $\{JR_{i}J\}_{i=1}^{N}$. Thus, we can recover either the original molecule or its mirror image.

The most commonly used procedure for estimating the common line of a pair of images $P_{R_{i}}$ and $P_{R_{j}}$, is to calculate their Fourier transforms and then find a maximally correlated pair of central lines between the transformed images, see e.g. [12]. As the images obtained by cryo-EM are contaminated with high levels of noise, in practice the probability of correctly detecting their common lines, and subsequently their rotations using the angular reconstitution method, is low. In [11], Shkolnisky and Singer describe an algorithm for estimating the rotations which achieves robustness to noise by employing an approach known as ’synchronization’. In this approach, the rotations $\{R_{i}\}_{i=1}^{N}$ are estimated using all the relative rotations $\{R_{i}^{T}R_{j}\}_{i<j\in[N]}$ together at once. The authors use the set $\{R_{i}^{T}R_{j}\}_{i<j\in[N]}$ to construct a $3N\times 3N$ block matrix $M$ known as a ’synchronization matrix’, whose $(i,j)^{th}$ block $M_{ij}$ of size $3\times 3$ is given by $R_{i}^{T}R_{j}$, that is

[TABLE]

Defining the matrix $U=\left(R_{1},\ldots,R_{N}\right)$, we see that

[TABLE]

and thus we can obtain $U=(R_{1},\ldots,R_{N})$ by factoring $M$ using SVD.

A method for estimating the set of relative rotations $\{R_{i}^{T}R_{j}\}_{i<j\in[N]}$ in a non-sequential manner is given in [12], and takes advantage of the following observation: the relative rotation of $P_{R_{i}},P_{R_{j}}$ can be estimated from the common lines between these images and any of the $N-2$ images $P_{R_{k}}$, where $k\neq i,j$. The authors of [12] show how to obtain a robust estimate of the relative rotation of each pair $P_{R_{i}}$ and $P_{R_{j}}$, by taking a majority vote over all these $N-2$ estimates. Note that in order to construct the synchronization matrix $M$, one has to obtain a hand-consistent set of relative rotations, i.e. either the set $\{R_{i}^{T}R_{j}\}_{i<j\in[N]}$ or the set $\{JR_{i}^{T}R_{j}J\}_{i<j\in[N]}$. However, if the estimates $R_{i}^{T}R_{j}$ or $JR_{i}^{T}R_{j}J$ are obtained independently for each pair of images, this cannot be guaranteed. A solution to this issue is given in [8].

A direct application of any of the common lines based methods described above to a symmetric molecule encounters substantial difficulties stemming from the ambiguity described by (2.7). Suppose that given a pair of images $P_{R_{i}}$ and $P_{R_{j}}$ of a $D_{2}$-symmetric molecule, we wish to estimate the relative rotation $R_{i}^{T}R_{j}$ by applying the angular reconstitution method. By (2.7), we can detect four different common lines between the images, corresponding to four different pairs of projection planes, and we have no way of knowing which common line corresponds to which pair of planes. Thus, combining the common line of $P_{R_{i}}$, $P_{R_{j}}$ together with common lines with $P_{R_{k}}$, gives rise to $4^{3}$ different possible combinations of common line triplets, many of which do not submit an intersection of three planes. For instance, we can erroneously consider a combination of the common lines between the pairs $\{P_{R_{i}},P_{R_{j}}\}$, $\{P_{R_{i}},P_{R_{k}}\}$, and a third pair $\{P_{R_{j}},P_{g_{2}R_{k}}\}$, from which we cannot establish the relative rotations of $P_{R_{i}},P_{R_{j}}$ and $P_{R_{k}}$, since we are not considering the correct common lines triplet between the projection planes of these images.

In [12], the authors derive a simple condition by which one can determine whether a triplet of common lines can be realized as the intersection of three central planes. Still, even common line triplets which do satisfy this condition can generate any of the rotations $\{R_{i}^{T}R_{j},R_{i}^{T}g_{2}R_{j},R_{i}^{T}g_{3}R_{j},R_{i}^{T}g_{4}R_{j}\}$, between which we cannot distinguish. Thus, to construct the synchronization matrix $M$ in (3.2), one would have to devise a way to obtain a set of estimates $\{\widetilde{R}_{i}^{T}\widetilde{R}_{j}\}_{i<j\in[N]}$, in which for each $i\in[N]$ all the relative rotations $\widetilde{R}_{i}^{T}\widetilde{R}_{j}$ for $j\neq i$ collectively ’agree’ on the identity of $\widetilde{R}_{i}\in\{R_{i},g_{2}R_{i},g_{3}R_{i},g_{4}R_{i}\}$.

A further difficulty stems from the fact that pairs of central lines which are adjacent to each common line between a pair of images are also highly correlated. Thus, attempting to estimate 4 common lines between a pair of noisy images $\hat{P}_{R_{i}}$ and $\hat{P}_{R_{j}}$ by simply trying to detect the 4 best correlated central lines is most likely to fail (see [9] for a detailed discussion).

In Section 4 we present a different approach for estimating all four common lines and respective relative rotations of a pair of projection images of a $D_{2}$-symmetric molecule, inspired by maximum likelihood methods. In Section 5 we outline an algorithm for extracting the rotations $R_{i}$ of $\eqref{Intro:ProjEq}$ from the relative rotations estimated by the procedure described in Section 4. This approach encounters 3 major obstacles which are resolved in Sections 6, 7 and 8. In Section 9 we demonstrate the applicability of our method to experimental cryo-EM data. Finally, in Section 10 we summarize and discuss future work.

4 Relative rotations estimation

In this section we present a method for estimating the set of relative rotations $\{R_{i}^{T}g_{m}R_{j}\}_{m=1}^{4}$ for a pair of images $P_{R_{i}}$ and $P_{R_{j}}$ of a $D_{2}$-symmetric molecule. Similarly to the maximal correlations approach described in the previous section, we begin by computing the 2D Fourier transform of each image $P_{R_{i}}$. By (2.7), in the noiseless case, each pair of transformed images $\hat{P}_{R_{i}}$ and $\hat{P}_{R_{j}}$ has exactly four pairs of perfectly correlated central lines. Thus, in principle, we can detect these common lines by computing correlations between pairs of central lines in $\hat{P}_{R_{i}}$ and $\hat{P}_{R_{j}}$, and choosing the four maximally correlated pairs. However, as was explained in the previous section, this approach encounters several substantial difficulties. We now describe a different approach inspired by maximum likelihood methods.

Consider the set

[TABLE]

of quadruplets of relative rotations generated from all pairs of rotations $Q_{l},Q_{r}\in SO(3)$ with non-coinciding beaming directions, and let us denote the members of $\mathcal{D}_{c}$ by $Q_{lr}=\{Q_{l}^{T}g_{m}Q_{r}\}_{m=1}^{4}$. Given a pair of images $P_{R_{i}}$ and $P_{R_{j}}$, we now show how one can use the common lines of the images to assign a score $\pi_{ij}(Q_{l},Q_{r})$ to each element $Q_{lr}\in\mathcal{D}_{c}$, which indicates how well it approximates the quadruplet $\{R_{i}^{T}g_{m}R_{j}\}_{m=1}^{4}$. Since $\{R_{i}^{T}g_{m}R_{j}\}_{m=1}^{4}\in\mathcal{D}_{c}$, it can be detected by searching over $D_{c}$ for a candidate $Q_{lr}$ with the best score $\pi_{ij}(Q_{l},Q_{r})$. We henceforth refer to $D_{c}$ as the relative rotations search space.

First, to relate each candidate $Q_{lr}\in\mathcal{D}_{c}$ to the common lines of $P_{R_{i}}$ and $P_{R_{j}}$, let us compute the vectors

[TABLE]

analogously to (2.5). If $Q_{lr}=\{R_{i}^{T}g_{m}R_{j}\}_{m=1}^{4}$, then the set $\{\tilde{q}_{lr}^{m}\}_{m=1}^{4}$ corresponds to the direction vectors of the common lines of $P_{R_{i}}$ and $P_{R_{j}}$. We subsequently refer to $\{\tilde{q}_{lr}^{m}\}_{m=1}^{4}$ as the set of common lines directions of the quadruplet $Q_{lr}$ (corresponding to the pair $Q_{l},Q_{r}\in SO(3)$). Next, for each candidate $Q_{lr}$, we use the set $\{\tilde{q}_{lr}^{m}\}_{m=1}^{4}$ to compute the coordinate vectors

[TABLE]

for $m=1,2,3,4$, analogously to (2.6). If $Q_{lr}=\{R_{i}^{T}g_{m}R_{j}\}_{m=1}^{4}$, the coordinates in (4.3) correspond to the local coordinates of the common lines of $P_{R_{i}}$ and $P_{R_{j}}$ on the respective projection planes of the images. Denote by

[TABLE]

the half line (known as a Fourier ray) in the direction which forms an angle $\theta$ with the $x$-axis of the transformed image $\hat{P}_{R_{n}}$. We then compute the normalized cross correlations

[TABLE]

of each pair of rays given by the direction vectors in (4.3). We use rays instead of lines, since the correlation $\rho_{ij}(\theta,\varphi)$ between each pair of rays $\nu_{i,\theta}(\xi)$ and $\nu_{j,\varphi}(\xi)$ is identical to the correlation value $\rho_{ij}(\theta+\pi,\varphi+\pi)$ of their anti-podal rays. We then assign to each quadruplet $Q_{lr}$ the score

[TABLE]

By (2.7), if $Q_{lr}=\{R_{i}^{T}g_{m}R_{j}\}_{m=1}^{4}$ for some $Q_{l},Q_{r}\in SO(3)$, then $\pi_{ij}(Q_{l},Q_{r})=1$. Thus, the quadruplet $Q_{lr}$ with $\pi_{ij}(Q_{l},Q_{r})=1$ is declared as $\{R_{i}^{T}g_{m}R_{j}\}_{m=1}^{4}$.

We remark that since in practice we use a discretization of $\mathcal{D}_{c}$ to estimate the relative rotations of pairs of noisy images (see Section 9.1 for details), $\pi_{ij}$ is never exactly 1, and so we simply choose a candidate $Q_{lr}$ which maximizes $\pi_{ij}(Q_{l},Q_{r})$ as an approximation for $\{R_{i}^{T}g_{m}R_{j}\}_{m=1}^{4}$.

However, we can obtain more robust estimates to $\{R_{i}^{T}g_{m}R_{j}\}_{m=1}^{4}$ by also combining self common lines into the score (4.5). As was explained in Section 2, each image $P_{R_{i}}$ of a $D_{2}$-symmetric molecule has three self common lines, given by (2.9), which are the intersections of the projection plane of $P_{R_{i}}$ with the projection planes of $P_{g_{2}R_{i}}$, $P_{g_{3}R_{i}}$ and $P_{g_{3}R_{i}}$. We now show how to adjust the score $\pi_{ij}(Q_{l},Q_{r})$ of each candidate $Q_{lr}\in\mathcal{D}_{c}$ to account for the self common lines of each image in the pair $P_{R_{i}}$ and $P_{R_{j}}$.

For each candidate $Q_{lr}\in\mathcal{D}_{c}$, we first compute the vectors

[TABLE]

analogously to (2.8). If $Q_{lr}=\{R_{i}^{T}g_{m}R_{j}\}_{m=1}^{4}$, then the sets $\{\tilde{q}_{ll}^{m}\}_{m=2}^{4}$ and $\{\tilde{q}_{rr}^{m}\}_{m=2}^{4}$ correspond to the directions of the self common lines of the images $P_{R_{i}}$ and $P_{R_{j}}$, respectively (see (2.8)). Next, we use $\{\tilde{q}_{ll}^{m}\}_{m=2}^{4}$ and $\{\tilde{q}_{rr}^{m}\}_{m=2}^{4}$ of (4.6) to compute the coordinates

[TABLE]

for $m=2,3,4$, analogously to (2.10). If $Q_{lr}=\{R_{i}^{T}g_{m}R_{j}\}_{m=1}^{4}$, then the coordinates in (4.7) correspond to the local coordinates of the self common lines of $P_{R_{i}}$ and $P_{R_{j}}$, on their respective projection planes (see (2.10)). We then compute the set of normalized autocorrelations

[TABLE]

for $m\in\{2,3,4\}$, and if $Q_{lr}=\{R_{i}^{T}g_{m}R_{j}\}_{m=1}^{4}$ for some $Q_{l},Q_{r}\in SO(3)$, then by (2.9) we have $\prod_{m=2}^{4}\rho_{ii}(\tilde{\alpha}_{ll}^{1m},\tilde{\alpha}_{ll}^{m1})\rho_{jj}(\tilde{\alpha}_{rr}^{1m},\tilde{\alpha}_{rr}^{m1})=1$. We therefore redefine the score $\pi_{ij}(Q_{l},Q_{r})$ in (4.5) to be

[TABLE]

For each $i<j\in[N]$, we then set

[TABLE]

to be $\{R_{i}^{T}g_{m}R_{j}\}_{m=1}^{4}$.

We remark, that as was explained in Section 1 and Section 3 (see (1.6)), the images $\hat{P}_{R_{n}}$ and $\hat{P}_{JR_{n}J}$ are identical. Thus, the self common lines of each image $\hat{P}_{R_{n}}$ are identical to the self common lines of $\hat{P}_{JR_{n}J}$, and the common lines between each pair of images $\hat{P}_{R_{i}}$ and $\hat{P}_{R_{j}}$ are identical to the common lines between $\hat{P}_{JR_{i}J}$ and $\hat{P}_{JR_{j}J}$. Hence, for each $i<j\in[N]$, the set $\{Q_{l}^{T}g_{m}Q_{r}\}_{m=1}^{4}\in\mathcal{D}_{c}$ which maximizes the score $\pi_{ij}$ has the same score as $\{JQ_{l}^{T}g_{m}Q_{r}J\}_{m=1}^{4}\in\mathcal{D}_{c}$. Thus, in (4.10), for each $i<j\in[N]$, we either estimate $\{R_{i}^{T}g_{m}R_{j}\}_{m=1}^{4}$ or $\{JR_{i}^{T}g_{m}R_{j}J\}_{m=1}^{4}$, independently from other pairs of $i$ and $j$.

The procedure for estimating of the sets of relative rotations for each $i<j\in[N]$ is summarized in Algorithm 1.

5 Estimating the rotation matrices

In the previous section, we have shown how to estimate for each $i<j\in[N]$ either $\{R_{i}^{T}g_{m}R_{j}\}_{m=1}^{4}$ or $\{JR_{i}^{T}g_{m}R_{j}J\}_{m=1}^{4}$. In Section 6 below, we will show how to resolve the handedness ambiguity. Thus, in this section we will assume w.l.o.g that we have the sets $\{R_{i}^{T}g_{m}R_{j}\}_{m=1}^{4}$ for all $i<j\in[N]$, and outline how to recover the rotations $R_{i}$ in (1.2) from these sets.

To recover the matrices $R_{i}$ row by row, we use the following observation. For each $m\in\{1,2,3\}$, denote by $I_{m}$ the $3\!\times\!3$ diagonal matrix

[TABLE]

and note that

[TABLE]

where $g_{l}$ were defined in (1.3). Thus, for any pair of matrices $R_{i}$ and $R_{j}$ we have that

[TABLE]

for $m=1,2,3$, where $v_{i}^{m}$ and $v_{j}^{m}$ are the $m^{\text{th}}$ rows of the matrices $R_{i}$ and $R_{j}$, respectively. We can also compute the matrices $(v_{i}^{m})^{T}v_{i}^{m}$ for $m\in\{1,2,3\}$ and $i\in[N]$, by noting that since $v_{j}^{m}$ are rows of orthogonal matrices we have

[TABLE]

Since in practice the matrices $(v_{i}^{m})^{T}v_{j}^{m}$ are estimated from noisy images, we get a more robust estimate for $(v_{i}^{m})^{T}v_{i}^{m}$ by using all $j$, that is, by setting

[TABLE]

for $m\in\{1,2,3\}$ and $i\in[N]$.

Next, for each $m\in\{1,2,3\}$, we construct the $3N\times 3N$ matrix $H_{m}$ whose $(i,j)^{th}$ $3\times 3$ block is given by the rank 1 matrix $(v_{i}^{m})^{T}v_{j}^{m}$, and note that

[TABLE]

That is, $H_{m}$, $m=1,2,3$, are rank 1 matrices. We can now factorize each matrix $H_{m}$ using SVD, to obtain either the vector $v_{m}$ or $-v_{m}$, hence retrieving either the set of rows $\{v_{i}^{m}\}_{i=1}^{N}$ or $\{-v_{i}^{m}\}_{i=1}^{N}$, for each $m\in\{1,2,3\}$. Then, we can use these sets of rows to assemble the matrices $\{OR_{i}\}_{i=1}^{N}$ row by row, where $O\in O(3)$ is a diagonal matrix with $\pm 1$ on its diagonal. If $\det{O}=-1$, we simply multiply all $OR_{i}$ by $-1$, and thus, we can assume w.l.o.g that $O$ is a rotation. The matrix $O$ is an inherent degree of freedom, since we can always “rotate the world” by any orthogonal matrix.

Unfortunately, the approach just described is not directly applicable, as we now explain. Recall from Section 3, that though we can recover the set of relative rotation matrices $\{R_{i}^{T}g_{m}R_{j}\}_{m=1}^{4}$ from the common lines of $P_{R_{i}}$ and $P_{R_{j}}$, we have no way of knowing for each $m\in\{1,2,3,4\}$ which of the recovered matrices in the latter set is $R_{i}^{T}g_{m}R_{j}$. This implies, that for each $i<j\in[N]$, we can only obtain a permutation $(R_{i}^{T}g_{\tau_{ij}(m)}R_{j})_{m=1}^{4}$ of the 4-tuple $(R_{i}^{T}g_{m}R_{j})_{m=1}^{4}$ where $\tau_{ij}\in S_{4}$ is some unknown permutation of $(1,2,3,4)$. In (5.3), we computed the 3-tuples $((v_{i}^{m})^{T}v_{j}^{m})_{m=1}^{3}$ by summing the first element of the 4-tuple $(R_{i}^{T}g_{m}R_{j})_{m=1}^{4}$ with the rest of its elements. Suppose for example, that we have the 4-tuple $(R_{i}^{T}g_{2}R_{j},R_{i}^{T}g_{3}R_{j},R_{i}^{T}g_{1}R_{j},R_{i}^{T}g_{4}R_{j})$ for a given a pair of images $P_{R_{i}}$ and $P_{R_{j}}$. One can easily verify by direct calculation that

[TABLE]

where $S_{3}$ is the group of all permutations of a 3-tuple. Now, observe that by (5.2) and (5.5) we have

[TABLE]

This implies that the summation in (5.3) of a permutation of the 4-tuple $(R_{i}^{T}g_{m}R_{j})_{m=1}^{4}$ results in a permutation of the respective 3-tuple $((v_{i}^{m})^{T}v_{j}^{m})_{m=1}^{3}$ of rank 1 matrices, where some of the matrices have a spurious $-1$ factor. The following proposition, the proof of which is given Appendix A.1, summarizes the effect of the aforementioned summation on a general permutation $(R_{i}^{T}g_{\tau(m)}R_{j})_{m=1}^{4}$ of a 4-tuple $(R_{i}^{T}g_{m}R_{j})_{m=1}^{4}$, that is, on the order and signs of the respective 3-tuple $((v_{i}^{m})^{T}v_{j}^{m})_{m=1}^{3}$.

Proposition 5.1**.**

Let $(R_{i}^{T}g_{\tau(m)}R_{j})_{m=1}^{4}$ for some $\tau\in S_{4}$ be a permutation of the $4$-tuple $(R_{i}^{T}g_{m}R_{j})_{m=1}^{4}$.

If $\tau(1)=1$, then the corresponding 3-tuple of rank 1 matrices is given by $((v_{i}^{m_{r}})^{T}v_{j}^{m_{r}})_{r=1}^{3}$, where $(m_{1},m_{2},m_{3})=(\tau(2)-1,\tau(3)-1,\tau(4)-1)$. 2. 2.

If $\tau(m)=1$ for $m>1$, then the corresponding 3-tuple of rank 1 matrices is given by

[TABLE]

Proposition 5.1 implies that we can only obtain the ordered triplets

[TABLE]

where the permutations $\sigma_{ij}\in S_{3}$ are unknown, and each of the matrices $(v_{i}^{\sigma_{ij}(m)})^{T}v_{j}^{\sigma_{ij}(m)}$ is multiplied by $\pm 1$, which is also unknown and depends on $\tau_{ij}$. Thus, we cannot construct the matrices $H_{m}$ in (5.4) using the triplets in (5.6) directly. We will show how to construct $H_{m}$ using the triplets (5.6) in Sections 7 and 8.

The following three sections are organized as follows. In Section 6, we show a method for handedness synchronization for $D_{2}$-symmetric molecules, which is adapted from a method proposed in [8] for non-symmetric molecules. In Section 7, we show how to partition the 3-tuples in (5.6) into three sets of the form $\{s_{ij}^{m}(v_{i}^{m})^{T}v_{j}^{m}\}_{i<j\in[N]}$ for $m=1,2,3$ and some unknown signs $s_{ij}\in\{-1,1\}$. Then, in Section 8 we show how to correct the signs $s_{ij}^{m}$ so that we can construct the $3N\times 3N$ matrices

[TABLE]

where $s_{i}^{m}\in\{-1,1\}$ for $i\in[N]$ and $m\in\{1,2,3\}$. We can then factor the matrices $\widetilde{H}_{m}$ using SVD to obtain the vectors $v_{m}$ (of length $3N$), which give us the sets of rows $\{s_{i}^{m}v_{i}^{m}\}_{i=1}^{N}$ for $m\in\{1,2,3\}$, and assemble the matrices

[TABLE]

Since for each $i\in[N]$ either $D_{i}$ or $-D_{i}$ is in $\{g_{m}\}_{m=1}^{4}$, that is, $\hat{R}_{i}=\pm g_{m}R_{i}$ for some $m\in\{1,2,3,4\}$, we can compute the matrices

[TABLE]

by replacing all matrices $\hat{R}_{i}$ which have $\det(\hat{R}_{i})=-1$ with $-\hat{R}_{i}$. The resulting set of matrices $\{\widetilde{R}_{i}\}_{i=1}^{N}$ satisfies $\widetilde{R}_{i}\in\{g_{m}R_{i}\}_{m=1}^{4}$ for all $i\in[N]$, and are therefore a solution for the orientation assignment problem which was stated at the end of Section 1.

6 Handedness synchronization

Following the discussion in the previous section, we now assume we have obtained a set of 4-tuples

[TABLE]

for some unknown $\delta_{ij}\in\{0,1\}$, by applying Algorithm 1. We now explain how to extract one of the hand-consistent sets

[TABLE]

from the set in (6.1).

For all $i<j\in[N]$, we denote by

[TABLE]

a 4-tuple of relative rotations consistent with a pair of images $P_{R_{i}}$ and $P_{R_{j}}$ of a $D_{2}$-symmetric molecule. We also denote by $R_{ij}^{m}$, the $m^{th}$ relative rotation in the 4-tuple $R_{ij}$ for $m=1,2,3,4$, and by $JR_{ij}J$ the set $\{JR_{ij}^{m}J\}_{m=1}^{4}$. We now show how the set in (6.1) can be partitioned into two disjoint sets

[TABLE]

Once we have the partition in (6.4), we can compute the set of 4-tuples $\widetilde{C}_{1}=\{(JR_{ij}^{m}J)_{m=1}^{4}|R_{ij}\in C_{1}\}$. Then, one of the hand-consistent sets in (6.2) is given by $C_{0}\cup\widetilde{C}_{1}$.

The partition in (6.4) is derived by the following procedure. First, we construct a graph $\Sigma$ with vertices corresponding to the estimates $R_{ij}$ in (6.3), and with edges that encode which pairs of estimates $R_{ij}$ and $R_{kl}$ are in the same set in (6.4), and which aren’t (as will be explained shortly). Then, we derive the partition in (6.4) from the eigenvector of the leading eigenvalue of the adjacency matrix of $\Sigma$. The procedure we present is an adaptation of an algorithm that was derived in [8] for non-symmetric molecules.

We next state a proposition, the proof of which is given in Appendix A.2, which allows us to determine which estimates $R_{ij}$, $i<j\in[N]$, belong to the same set in (6.4). Following the approach in [8], we look at triplets of estimates of the form $R_{ij},R_{jk},R_{ki}$ for all triplets $i<j<k\in[N]$, and determine which members of each such triplet are in the same set of (6.4) and which aren’t.

Proposition 6.1**.**

For any $i<j<k\in[N]$ consider the triplet of estimates

[TABLE]

that is, $R_{ij},R_{jk}$ and $R_{ki}$ are all in the same set of (6.4). Then, exactly 16 of the $4^{3}$ matrix products in the set

[TABLE]

satisfy

[TABLE]

Now, consider a triplet of estimates

[TABLE]

and note that since each estimate is either in the set $C_{0}$ or in the set $C_{1}$ of (6.4), it must be that either all estimates are in the same set, or two estimates are in one set and the third estimate is in the other. We define the “set configuration” of a triplet $(R_{ij},R_{jk},R_{ki})$ by the row vector

[TABLE]

and denote

[TABLE]

Note that if, for example, $R_{ij}$ is in one set of (6.4) and $R_{jk}$ and $R_{ki}$ are in another, then we have that $JR_{ij}J,R_{jk}$ and $R_{ki}$ are all in the same set. We remark that we have found experimentally that whenever three estimates $R_{ij}$,$R_{jk}$ and $R_{ki}$ are not in the same set of (6.4), then all the products in (6.7) are far from $I$ in norm. Thus, Proposition 6.1 suggests that we can find the set configuration of a triplet of estimates by the following procedure. First, we compute the four sets of $4^{3}$ norms

[TABLE]

where $\|\cdot\|_{F}$ is the Frobenius norm. Next, we sort the norms in each set $\mathcal{N}_{ijk}^{0},\mathcal{N}_{ijk}^{1},\mathcal{N}_{ijk}^{2}$ and $\mathcal{N}_{ijk}^{3}$ in (6.11) in ascending order, and denote the resulting ascending sequences by $\mathcal{S}_{ijk}^{0},\mathcal{S}_{ijk}^{1},\mathcal{S}_{ijk}^{2}$ and $\mathcal{S}_{ijk}^{3}$, respectively. Finally, we compute the scores

[TABLE]

By Proposition 6.1, there are exactly 16 norms with value 0 in the set $\mathcal{N}_{ijk}^{p}$ of (6.11) which corresponds to the correct set configuration $d_{ijk}$ of a triplet $(R_{ij},R_{jk},R_{ki})$ (see (6.9)). Thus, we set $d_{ijk}=c_{p}$ for $p\in\{0,1,2,3\}$ such that $\hat{\mathcal{S}}_{ijk}^{p}$ is the minimal score in (6.12).

Once we have computed $d_{ijk}$ for all $i<j<k\in[N]$, we construct a graph $\Sigma$ whose vertices correspond to the estimates $R_{ij}$ in (6.3), and whose edges are defined by the $\binom{N}{2}\times\binom{N}{2}$ adjacency matrix (which we also denote by $\Sigma$)

[TABLE]

Finally, we compute the eigenvector $u_{s}$ which corresponds to the leading eigenvalue of the matrix $\Sigma$. In [8], it is shown that the leading eigenvalue of $\Sigma$ is simple, and that $u_{s}$ is of the form $\{-1,1\}^{\binom{N}{2}}$ (up to normalization), where the sign of each entry encodes the set membership in (6.4) of each estimate $R_{ij}$. The procedure for handedness synchronization for $D_{2}$-symmetric molecules is summarized in Algorithm 2.

7 Rotations’ rows synchronization

At this point, in light of Sections 5 and 6, we assume that we have obtained one of the hand-consistent sets of 4-tuples in (6.2). Let us assume without loss of generality that we have the set $\{(R_{i}^{T}g_{\tau_{ij}(m)}R_{j})_{m=1}^{4}\}_{i<j\in[N]}$. As was explained in Section 5, for each $i<j\in[N]$, we now form a 3-tuple of matrices by summing the first element of $(R_{i}^{T}g_{\sigma_{ij}(m)}R_{j})_{m=1}^{4}$ with each of the rest of its elements. By Proposition 5.1, this results in a set of triplets

[TABLE]

which was defined in (5.6), where $\sigma_{ij}$ are unknown and the $\pm 1$ signs are also unknown. In this section, we will show how to partition this set of triplets into three disjoint sets

[TABLE]

where $s_{ij}^{m}$ are the (unknown) signs of $(v_{i}^{m})^{T}v_{j}^{m}$. That is, for each $m\in\{1,2,3\}$, the set $C_{m}$ contains all outer products between the $m^{th}$ rows of the rotation matrices $R_{i}$ and $R_{j}$ for $i<j\in[N]$, up to sign. This partition will be obtained by casting it as a graph partitioning problem.

In Section 7.1, we show how to encode the partition in (7.1) as a graph in which each vertex corresponds to one of the matrices in (5.6). In Section 7.2, we construct the adjacency matrix of the graph, and in Section 7.3, we show how to extract the partition in (7.1) from the leading eigenvectors of the graphs’ adjacency matrix.

7.1 Graph partitioning formulation

In what follows, we denote the $3\times 3$ matrices in (7.1) by

[TABLE]

We now construct a weighted graph $\Omega=(V,E)$ from which the partition in (7.1) can be inferred. Each vertex in $V$ corresponds to one of the matrices in (7.2) (henceforth, we shall refer to both the matrix $s_{ij}^{m}(v_{i}^{m})^{T}v_{j}^{m}$ and its corresponding vertex in $V$ using the notation $v_{ij}^{m}$). Thus, we have that (see (7.2))

[TABLE]

We define the set of weighted edges $E$ of the graph $\Omega$ by its $3\binom{N}{2}\times 3\binom{N}{2}$ adjacency matrix, which we also denote by $\Omega$, as follows

[TABLE]

That is, we only connect by an edge vertices which have exactly one index in common. We give this edge a weight $+1$ if its incident vertices are in the same set $C_{m}$ of (7.3), and weight $-1$ otherwise (see Fig. 1).

Note that the weights on the edges $E$ of $\Omega$ induce a partition of the vertex set $V$ into the sets of (7.3), by grouping together vertices which are connected by edges with a weight of $+1$. Recovering the partition in (7.3) corresponds to coloring the vertices $V$ of $\Omega$ with 3 colors, say, red, green, and blue, where the vertices of $C_{1}$ are colored red, of $C_{2}$ green and of $C_{3}$ blue.

In what follows, we show that the partition of $V$ to the sets $C_{m}$ in (7.3) can be derived from eigenvectors of the matrix $\Omega$.

Definition 7.1**.**

Let the eigenvalues of an $n\times n$ matrix A be $\lambda_{1}>\lambda_{2}>\ldots>\lambda_{r}$, with their respective multiplicities given by $n_{1},\ldots,n_{r}$. We denote the spectrum of A by $\Lambda(A)$, and write

[TABLE]

Theorem 7.2**.**

The spectrum of the matrix $\Omega$ is given by

[TABLE]

The proof of Theorem 7.2 is given in Appendix A.3.

Definition 7.3**.**

Define $\alpha=(2{\binom{N}{2}})^{-\frac{1}{2}}$, $\beta=(6{\binom{N}{2}})^{-\frac{1}{2}}$, and define the pair of vectors $u_{\alpha},u_{\beta}\in\mathbb{R}^{3{\binom{N}{2}}}$ by

[TABLE]

where for any $w\in\mathbb{R}^{3{\binom{N}{2}}}$ we denote by $w(v_{ij}^{m})$ the entry of $w$ with the same index as the row of $\Omega$ which corresponds to the vertex $v_{ij}^{m}$.

Throughout Section 7, for any column vector $w\in\mathbb{R}^{3{\binom{N}{2}}}$, we denote

[TABLE]

That is, $(w)_{ij}\in\mathbb{R}^{3}$ is the column vector that corresponds to the entries of the triplet $(v_{ij}^{\sigma_{ij}(1)},v_{ij}^{\sigma_{ij}(2)},v_{ij}^{\sigma_{ij}(3)})$ in $w$.

Proposition 7.4**.**

The vectors $u_{\alpha}$ and $u_{\beta}$ in Definition 7.3 are orthogonal eigenvectors of $\Omega$, corresponding to the eigenvalue $\mu_{c}=4(N-2)$.

We prove Proposition 7.4 in Appendix A.4. The immediate consequence of Theorem 7.2 and Proposition 7.4 is the following corollary.

Corollary 7.5**.**

The eigenspace of $\mu_{c}$ is spanned by $u_{\alpha}$ and $u_{\beta}$.

Note that $u_{\alpha}$ is a unit vector which exactly encodes the partition in (7.3), where the entries $+\alpha,0$ and $-\alpha$ encode the color of each vertex $v_{ij}^{m}$. The unit vector $u_{\beta}$ is ’color blind’ in the sense that it can only distinguish between 2 colors. Obviously, in any 3-coloring of a graph we can always permute the colors, e.g., switch the color of all red vertices to green, green vertices to blue, and blue vertices to red. This is manifested in the following proposition.

Proposition 7.6**.**

Define the ’3-color’ and ’2-color’ vectors by

[TABLE]

respectively. For any $\sigma\in S_{3}$, we define the vectors $u_{\alpha}^{\sigma},u_{\beta}^{\sigma}\in\mathbb{R}^{3{\binom{N}{2}}}$ by

[TABLE]

where $u_{\alpha}^{\sigma}(v_{ij}^{m})$ and $u_{\beta}^{\sigma}(v_{ij}^{m})$ are the entries of $u_{\alpha}^{\sigma}$ and $u_{\beta}^{\sigma}$ with the same index as the row of $\Omega$ which corresponds to the vertex $v_{ij}^{m}$. Then, $u_{\alpha}^{\sigma}$ and $u_{\beta}^{\sigma}$ are orthogonal eigenvectors of $\Omega$ of (7.4) in the eigenspace of $\mu_{c}=4(N-2)$.

Proof.

Observe that $u_{\alpha}^{\sigma}$ is obtained from $u_{\alpha}$ of Definition 7.3 by replacing all $+\alpha$ entries with $u_{3c}(\sigma(1))$, all 0 entries with $u_{3c}(\sigma(2))$, and all $-\alpha$ entries with $u_{3c}(\sigma(3))$. The vector $u_{\beta}^{\sigma}$ is obtained from $u_{\beta}$ in a similar manner. The proposition then follows by repeating the method of proof applied in Proposition 7.4 with $u_{\alpha}^{\sigma}$ and $u_{\beta}^{\sigma}$.

∎

Following Corollary 7.5, we recover $u_{\alpha}$ (up to a color permutation, i.e, one of the vectors $u_{\alpha}^{\sigma}$ of (7.8)) in the following manner. We begin by constructing the matrix $\Omega$. Then, we compute a pair of orthogonal eigenvectors $v_{a}$ and $v_{b}$ spanning the eigenspace of $\mu_{c}$ (the leading eigenvalue of $\Omega$). In general, each of these eigenvectors is an orthogonal linear combination of $u_{\alpha}$ and $u_{\beta}$, and thus, we cannot read the partition in (7.3) directly from any one of them. In Section 7.3, we show how to ’unmix’ $v_{a}$ and $v_{b}$ and retrieve $u_{\alpha}$. In practice, due to noise, we can only compute an approximation of $u_{\alpha}$, and thus, we never get the exact values $\alpha$, [math] and $-\alpha$. We explain how to deal with this issue in Section 7.3. In the following section, we show how to construct $\Omega$ of (7.4) using the set of matrices in (5.6).

7.2 Constructing $\Omega$

We now derive a procedure for constructing the matrix $\Omega$ of (7.4). For any two pairs of indices $i<j\in[N]$ and $k<l\in[N]$, we denote by

[TABLE]

the $3\times 3$ matrix given by the rows of $\Omega$ corresponding to the vertices $v_{ij}^{\sigma_{ij}(1)},v_{ij}^{\sigma_{ij}(2)}$ and $v_{ij}^{\sigma_{ij}(3)}$, and columns of $\Omega$ corresponding to the vertices $v_{kl}^{\sigma_{kl}(1)},v_{kl}^{\sigma_{kl}(2)}$ and $v_{kl}^{\sigma_{kl}(3)}$. By (7.4) and (7.9), we have

[TABLE]

where $0_{3\times 3}$ is the $3\times 3$ zero matrix. We will now show how to construct $\Omega$ block by block, by computing the blocks $\Omega_{(i,j),(k,l)}$ for which $|\{i,j\}\cap\{k,l\}|=1$.

The following lemma, the proof of which is given in Appendix A.5, characterizes the indices of the non-zero entries in $\Omega$ of (7.4).

Lemma 7.7**.**

Define

[TABLE]

and for $i<j\in[N]$ define

[TABLE]

Moreover, for $i<j<k\in[N]$, define

[TABLE]

Then, we have that

[TABLE]

By (7.10) and the second equality in (7.14) of Lemma 7.7, to construct $\Omega$, we only need to determine the blocks

[TABLE]

which we now show how to do.

For each triplet of indices $i<j<k\in[N]$, we consider the vertices of $\Omega$ corresponding to the pairs of indices $(i,j)$, $(j,k)$, and $(i,k)$, written in the rows of the table

[TABLE]

For each pair of vertices from different rows in (7.16), we need to determine whether this pair belongs to the same set $C_{m}$ of (7.3) or to different sets. This corresponds to choosing between an edge with a weight of $+1$ or $-1$ for each of these pairs in $\Omega$. We therefore show how to determine all edge weights between the vertices in (7.16) simultaneously. This procedure is then repeated for each triplet of indices $i<j<k\in[N]$.

First, observe that for any pair of matrices $v_{ij}^{m}$ and $v_{jk}^{r}$, given by (7.2), we have

[TABLE]

since the row vectors $v_{j}^{m}$ and $v_{j}^{r}$ are rows of the orthogonal matrix $R_{j}$. This suggests that for each pair of vertices $v_{ij}^{\sigma_{jk}(m)}$ and $v_{jk}^{\sigma_{ij}(r)}$ with unknown $\sigma_{ij},\sigma_{jk}\in S_{3}$, we can determine whether they belong to the same set of (7.3), by simply computing the norm of the product of the matrices which they represent.

Since in practice we work with noisy data, we next show how to get more robust estimates for the edge weights of $\Omega$, by leveraging the graph structure of $\Omega$ in conjunction with (7.17). Denote by

[TABLE]

the transposed matrices of (7.2). By (7.17) and (7.18), for each triplet of matrices $v_{ij}^{m},v_{jk}^{r}$ and $v_{ki}^{p}=(v_{ik}^{p})^{T}$, $i<j<k\in[N]$, we have

[TABLE]

and by (7.18), we also have that

[TABLE]

Note that a non-zero product of matrices in (7.19) corresponds to a triplet of vertices in $\Omega$ in the same set of (7.3) (see Fig. 1(b)). We now infer which vertices in (7.16) belong to the same set of (7.3), by constructing a function that vanishes for all vertex triplets for which (7.19) and (7.20) hold, namely, for all vertex triplets that belong to the same set of (7.3). Specifically, for each triplet of indices $i<j<k\in[N]$, we minimize the function (which will be explained shortly) $f_{ijk}:S_{3}\times S_{3}\to\mathbb{R}$ given by

[TABLE]

over all $\gamma=(\gamma_{1},\gamma_{2},\gamma_{3})$ and $\delta=(\delta_{1},\delta_{2},\delta_{3})$ in $S_{3}$, and all choices of sign $\pm 1$ between the 2 terms in each norm (since by (7.19), the sign of the right term in each norm is unknown), independently between the norms.

The rationale of minimizing (7.21) can be demonstrated in the following manner. Writing down the vertices as in (7.16), we seek to rearrange the vertices in the second and third rows so that after rearrangement, the vertices in each column are in the same set of (7.3) (see example in Fig. 2). Whenever we choose a pair of permutations $\gamma^{*},\delta^{*}\in S_{3}$ such that the vertices in each column in

[TABLE]

are in the same set of (7.3), by (7.19) and (7.20), we have that all the terms in the sum (7.21) equal zero. Otherwise, if there exists a column in (7.22) in which there is a pair of vertices in different sets of (7.3), then we have that (7.21) is strictly $>0$. For example, if in each column of (7.22), there is a pair vertices in different classes, then by (7.19), the left term inside each norm in (7.21) equals zero, while the right term in each norm equals $(v_{ij}^{\sigma_{ij}(m)})^{T}v_{ji}^{\sigma_{ij}(m)}$, and we get that $f_{ijk}=\sum_{m=1}^{3}\|(v_{ij}^{\sigma_{ij}(m)})^{T}v_{ji}^{\sigma_{ij}(m)}\|$.

Once we compute for each $i<j<k\in[N]$ a pair of permutations $\gamma^{*}$ and $\delta^{*}$ minimizing $f_{ijk}$ in (7.21), the matrix $\Omega$ is set block by block by computing all the blocks of (7.15) in the following manner. We first set the first three blocks of (7.15), that is, $\Omega_{(i,j)(j,k)}$, $\Omega_{(i,j)(i,k)}$ and $\Omega_{(j,k)(i,k)}$. Consider (7.22), which consists of all the vertices that are incident to the edges that constitute the aforementioned blocks (see (7.9)). The triplet of vertices $v_{ij}^{\sigma_{ij}(m)},v_{jk}^{\sigma_{jk}(\gamma^{*}(m))}$ and $v_{jk}^{\sigma_{jk}(\delta^{*}(m))}$, which are all in the same column of (7.22), are in the same set $C_{\sigma_{ij}(m)}$ of (7.3), $m=1,2,3$. Thus, we assign a weight +1 to the edges of $\Omega_{(i,j)(j,k)}$, $\Omega_{(i,j)(i,k)}$, and $\Omega_{(j,k)(i,k)}$, which correspond to the pairs of vertices $(v_{ij}^{\sigma_{ij}(m)},v_{jk}^{\sigma_{jk}(\gamma^{*}(m))})$, $(v_{ij}^{\sigma_{ij}(m)},v_{ik}^{\sigma_{ik}(\delta^{*}(m))})$, and $(v_{jk}^{\sigma_{jk}(\gamma^{*}(m))},v_{ik}^{\sigma_{ik}(\delta^{*}(m))})$, for $m=1,2,3$. All the edges of $\Omega_{(i,j)(j,k)}$, $\Omega_{(i,j)(i,k)}$, and $\Omega_{(j,k)(i,k)}$ which correspond to pairs of the form $(v_{ij}^{\sigma_{ij}(m)},v_{jk}^{\sigma_{jk}(\gamma^{*}(r))})$, $(v_{ij}^{\sigma_{ij}(m)},v_{ik}^{\sigma_{ik}(\delta^{*}(r))})$, and $(v_{jk}^{\sigma_{jk}(\gamma^{*}(m))},v_{ik}^{\sigma_{ik}(\delta^{*}(r))})$ where $m\neq r$, are assigned a weight $-1$ since they are in a different sets of (7.3). As for the last three blocks of (7.15), that is, $\Omega_{(j,k)(i,j)}$, $\Omega_{(i,k)(i,j)}$, and $\Omega_{(i,k)(j,k)}$, by (7.4) we have that $\Omega(v_{ij}^{m},v_{kl}^{r})=\Omega(v_{kl}^{r},v_{ij}^{m})$ for all $i<j\in[N]$, $k<l\in[N]$, and $m,r\in\{1,2,3\}$. Thus, by (7.9) we have

[TABLE]

Thus, for every $i<j<k\in[N]$ it holds that

[TABLE]

and thus, $\Omega_{(j,k)(i,j)}$, $\Omega_{(i,k)(i,j)}$ and $\Omega_{(i,k)(j,k)}$ can be set according to (7.24) after we compute $\Omega_{(i,j)(j,k)},\Omega_{(i,j)(i,k)}$ and $\Omega_{(j,k)(i,k)}$.

The procedure for constructing $\Omega$ of (7.4) is summarized in Algorithm 3. In the next section, we turn to the task of unmixing the eigenvectors corresponding to the maximal eigenvalue of $\Omega$, in order to extract $u_{\alpha}$ of Definition 7.3.

7.3 Unmixing the eigenvectors of $\mu_{c}$

By Corollary 7.5, the eigenspace of $\Omega$ of (7.4) corresponding the eigenvalue $\mu_{c}=4(N-2)$ is spanned by $u_{\alpha}$ and $u_{\beta}$ of Definition 7.3. However, any orthogonal linear combination of $u_{\alpha}$ and $u_{\beta}$ is also an eigenvector, and so we can only compute two orthogonal eigenvectors which are linear combinations of $u_{\alpha}$ and $u_{\beta}$. In this section, we show how to ’unmix’ these linear combinations to retrieve $u_{\alpha}$.

Suppose that we have computed a pair of orthogonal unit eigenvectors

[TABLE]

spanning the eigenspace of $\mu_{c}$. Since $v_{a}$ and $v_{b}$ are unit vectors, by Proposition 7.4 we have

[TABLE]

that is, $(a_{1},a_{2})^{T}$ is also a unit vector. Similarly, we have that $\|(b_{1},b_{2})^{T}\|=1$. Furthermore, again by Proposition 7.4, we have

[TABLE]

i.e., the coefficients vectors $(a_{1},a_{2})^{T}$ and $(b_{1},b_{2})^{T}$ are unit orthogonal vectors. Denoting

[TABLE]

by (7.26) and (7.27) there exists an angle $\varphi\in[0,2\pi)$ such that either

[TABLE]

and by (7.25) and (7.29), we have that

[TABLE]

However, it can be easily verified that

[TABLE]

and thus, (7.30) can be written as

[TABLE]

Equation (7.32) suggests that we can if we can recover $\varphi$, we can unmix $v_{a}$ and $v_{b}$ and recover either $u_{\alpha}$ or $-u_{\alpha}$. By Definition 7.3, $-u_{\alpha}$ is obtained from $u_{\alpha}$ by switching places between all $+\alpha$ and $-\alpha$ values in $u_{\alpha}$, and so both $u_{\alpha}$ and $-u_{\alpha}$ encode the same partition in (7.3) of the vertices of $\Omega$. We now show how to find $\varphi$.

For any angle $\theta\in[0,2\pi)$, we write

[TABLE]

In addition, using the notation introduced in (7.6), for any vector $w\in\mathbb{R}^{3{\binom{N}{2}}}$ and for all $i<j\in[N]$, we define

[TABLE]

and we define by $d_{ij}(w)$ the value of $(w)_{ij}$ whose magnitude is between $M_{ij}$ and $m_{ij}$. Then, we define the function $f_{c}:[0,2\pi)\to\mathbb{R}$ by

[TABLE]

The following proposition, the proof of which is given in Appendix A.6, states that the minimum of (7.35) over $\theta\in[0,2\pi)$ is obtained at $\theta$ for which $v_{a}^{\theta}=\pm u_{\alpha}$ and $v_{b}^{\theta}=u_{\beta}$ (up to a permutation of the vectors $u_{\alpha}$ and $u_{\beta}$ as in (7.8)).

Proposition 7.8**.**

Out of all orthogonal pairs of unit vectors in the eigenspace of $\mu_{c}$, the 12 pairs of vectors $\{(\pm u_{\alpha}^{\sigma},u_{\beta}^{\sigma})\;|\;\sigma\in S_{3}\}$ are the unique minimizers of $f_{c}$ in (7.35), up to normalization.

Thus, if we denote by $(v_{a}^{\theta_{*}},v_{b}^{\theta_{*}})$ a minimizer of (7.35), then we declare $u_{\alpha}$ to be $v_{a}^{\theta^{*}}$. The partition in (7.3) is then read from $u_{\alpha}$. In practice, due to noise in the input data, $v_{a}^{\theta^{*}}$ never exactly equals $u_{\alpha}$, and so after we compute $v_{a}^{\theta^{*}}$, we threshold its entries according to

[TABLE]

The estimation of $u_{\alpha}$ is summarized in Algorithm 4.

8 Signs synchronization

Assuming we have obtained the partition in (7.1) (see also (7.3)), our final task is to adjust the signs $s_{ij}^{m}$ in each set $C_{m}$ of (7.1), so that we can construct the rank 1 matrices $\widetilde{H}_{m}$ of (5.7). As was explained in Section 5, the matrices $\widetilde{H}_{m}$ can then be decomposed to retrieve all the rows of the matrices $R_{i}$ in (1.1), which can then be assembled from their constituent rows. Since all assertions we derive in this section apply identically and independently to each set $C_{m}$ in (7.1), throughout this section we refer to a single set of rank 1 matrices

[TABLE]

by dropping the superscript $m$ that indicates the set $C_{m}$ in (7.1) to which the matrices belong. Our goal is therefore to estimate $v_{1},\ldots,v_{N}$ up to an arbitrary sign. To that end, we will adjust the signs $s_{ij}$ in (8.1) so that the matrix $\widetilde{H}$ of size $3N\times 3N$ whose $(i,j)^{th}$ block of size $3\times 3$ is given by (8.1) has rank 1. The leading eigenvector of $\widetilde{H}$ will then give the vectors $v_{1},\ldots,v_{N}$ as required. We next describe the “signs adjustment” procedure and the construction of $\widetilde{H}$.

The task of adjusting the signs $s_{ij}$ in (8.1) consists of three steps. The first step of the signs adjustment procedure is computing the matrices $v_{i}^{T}v_{i}$ for all $i\in[N]$ by observing that

[TABLE]

since $s_{ij}\in\left\{-1,1\right\}$ and so $s_{ij}^{2}=1$. For notational convenience, for all $i\in[N]$ we write $s_{ii}v_{i}^{T}v_{i}$ instead $v_{i}^{T}v_{i}$, since the ‘sign’ $s_{ii}$ of $v_{i}^{T}v_{i}$ equals 1. In principle, (8.2) allows to compute $s_{ii}v_{i}^{T}v_{i}$ for each $i\in[N]$ by using a single matrix $s_{ij}v_{i}^{T}v_{j}$ for some arbitrarily chosen $j\in[N]\backslash\{i\}$. However, since in practice the input data is noisy, we obtain more robust estimates for $v_{i}^{T}v_{i}$ by computing the averages

[TABLE]

followed by computing the best rank 1 approximation of (8.3) using SVD.

While $s_{ij}$ in (8.1) are defined only for $i<j$, for notational convenience we define $s_{ij}=s_{ji}$ whenever $i>j$, and as explained above $s_{ii}=1$. Thus, $s_{ij}$ are defined for all $i,j\in[N]$. We next outline steps 2 and 3 of the signs adjustment procedure, before giving a detailed description of these steps. In step 2 of the procedure, we construct $N$ rank 1 matrices $H_{1}^{s},\ldots,H_{N}^{s}$, which admit the decompositions

[TABLE]

for unknown $s_{nj}\in\{-1,1\}$, $j\in[N]$. For each $n\in[N]$, the matrix $H_{n}^{s}$ in (8.4) is constructed block by block from the matrices in (8.1) and (8.3), in such a way that each $H_{n}^{s}$ is a rank 1 matrix which admits the decomposition (8.4). This construction relies on Proposition 8.1 stated below. Each of the matrices $H_{n}^{s}$ in (8.4) can then be decomposed to recover all the rows $v_{1},\ldots,v_{N}$ up to signs $s_{nj}$. Thus, in theory, we could construct only one of these matrices, say $H_{1}^{s}$, and recover $\{s_{11}v_{1},\ldots,s_{1N}v_{N}\}$, which is our goal in this section. However, in practice, the estimates $s_{ij}v_{i}^{T}v_{j}$ of (8.1) contain errors since they were estimated from noisy images. Moreover, for each particular $n\in[N]$, the set of estimates $\{s_{nj}v_{n}^{T}v_{j}\}_{j\in[N]\backslash\{n\}}$ used to construct the matrix $H_{n}^{s}$ critically depends on the common lines of the single noisy image $P_{R_{n}}$ in (1.2) with each of the images $P_{R_{j}}$ for $j\in[N]\backslash\{n\}$, which due to the noise in the input images, may be highly inaccurate, leading to large errors in $(s_{11}v_{1},\ldots,s_{1N}v_{N})$ of (8.4).

Thus, to use all available data in estimating $v_{1},\ldots,v_{N}$ of (8.4), we first decompose all matrices $H_{n}^{s}$ in (8.4), which results in $N$ independent estimates $\{s_{nj}v_{n}\}_{j=1}^{N}$ for each row $v_{n}$. We then execute the third step of our signs adjustment procedure in which we use all estimates $v_{n}^{s}$ together (see (8.4)) to obtain a set of signs $\widetilde{s}_{ij}\in\{-1,1\}$, which allows us to adjust the signs $s_{ij}$ of (8.1) by multiplying each matrix $s_{ij}v_{i}^{T}v_{j}$ by $\widetilde{s}_{ij}$, such that the matrix $\widetilde{H}$ of size $3N\times 3N$ whose $(i,j)^{th}$ $3\times 3$ block is $(\widetilde{s}_{ij}s_{ij})v_{i}^{T}v_{j}$ has rank 1. This latter procedure exploits all available data at once, improving the robustness of the estimation of the matrices $R_{i}$ of (1.2) to noisy input data. This matrix $\widetilde{H}$ admits the decomposition

[TABLE]

for some unknown signs $s_{n}\in\{-1,+1\}$. Recalling that we dropped the index $m\in\{1,2,3\}$ from $\widetilde{H}_{m}$ of (5.7), and constructed $\widetilde{H}$ from the set in (8.1), we see that in fact we can construct $\widetilde{H}_{m}$ of (5.7) using $\{s_{ij}^{m}v_{i}^{T}v_{j}\}_{i<j\in[N]}$ for each $m\in\{1,2,3\}$, as required. We can then decompose each $\widetilde{H}_{m}$, and recover the rotation matrices $R_{i}$ of (1.2), as was explained in Section 5, which is our task in this paper. We now complete the details of the signs adjustment procedure described above, i.e., the construction of the matrices of (8.4) and (8.5).

The following proposition, the proof of which is given in Appendix A.7, is the basis for the construction of the matrices of (8.4).

Proposition 8.1**.**

Let $H$ be a $3N\times 3N$ matrix whose $(i,j)^{th}$ block of size $3\times 3$ is given by $s_{ij}v_{i}^{T}v_{j}$ of (8.1) if $i<j$, by its transpose if $i>j$, and by $v_{i}^{T}v_{i}$ of (8.2) if $i=j$. Then, H is rank 1 iff for each $n\in[N]$

[TABLE]

where as noted above, for $i>j$ we define $s_{ij}=s_{ji}$. Furthermore, whenever (8.6) holds, we have

[TABLE]

Now, suppose we wish to construct $H_{1}^{s}$ of (8.4). Proposition 8.1, and in particular (8.6), suggest that we can construct $H_{1}^{s}$ using the set in (8.1), by applying the following sign correction procedure. Recall that $v_{i}$ are rows of $3\times 3$ orthogonal matrices, and thus $v_{i}v_{i}^{T}=1$ for all $i\in[N]$. For each pair $1<i<j\in[N]$, we compute the norm

[TABLE]

and replace the matrix $s_{ij}v_{i}^{T}v_{j}$ with $-s_{ij}v_{i}^{T}v_{j}$ if (8.8) is greater than zero. Let $\{\hat{s}_{ij}v_{i}^{T}v_{j}\}_{i<j\in[N]}$ be the resulting set of rank 1 matrices after this signs correction procedure. By construction, it holds that

[TABLE]

Now, let $H_{1}^{s}$ be the $3N\times 3N$ matrix whose $(i,j)^{th}$ block of size $3\times 3$ is given by $\hat{s}_{ij}v_{i}^{T}v_{j}$. By (8.9) and Proposition 8.1, $H_{1}^{s}$ admits the decomposition

[TABLE]

The matrices $H_{2}^{s},\ldots,H_{N}^{s}$ are obtained in a similar manner. In practice, (8.8) is never exactly zero due to errors stemming from noise, and thus we also compute $\|(s_{i1}v_{i}^{T}v_{1})(s_{1j}v_{1}^{T}v_{j})+s_{ij}v_{i}^{T}v_{j}\|_{F}$, and replace $s_{ij}v_{i}^{T}v_{j}$ with $-s_{ij}v_{j}^{T}v_{j}$, if

[TABLE]

At this point, we can factor each of the matrices $\{H_{n}^{s}\}_{n=1}^{N}$ (e.g., using SVD), and obtain the set of vectors

[TABLE]

of (8.4), where $s_{n}\in\{-1,+1\}$ are unknown. For each $i\in[N]$ we have

[TABLE]

that is, we have $N$ estimates $\{s_{n}s_{ni}v_{i}\}_{n\in[N]}$ for the row $v_{i}$, where each estimate has an unknown sign $s_{n}s_{ni}$. This concludes step 2 of the signs adjustment procedure outlined above.

For the third and final step of the signs adjustment procedure, we define

[TABLE]

Since $v_{1},\ldots,v_{N}$ are unit row vectors, by (8.11) and (8.12), we have

[TABLE]

for all $i\neq j\neq k\in[N]$. Thus, we can obtain all the products $\widetilde{s}_{ij}\widetilde{s}_{jk}$ in (8.13) by taking dot products of the vectors in (8.11).

Now, suppose we computed the set $\{\widetilde{s}_{ij}\}_{i<j\in[N]}$ from the values $\widetilde{s}_{ij}\widetilde{s}_{jk}$ in (8.13) (as will be explained shortly). Since $\widetilde{s}_{ij}=s_{i}s_{j}s_{ij}$, we have that

[TABLE]

Thus, we can multiply each matrix $s_{ij}(v_{i})^{T}v_{j}$ in (8.1) by $\widetilde{s}_{ij}$, and obtain the set of matrices $\{s_{i}s_{j}v_{i}^{T}v_{j}\}_{i<j\in[N]}$, and together with (8.3), we can construct the $3N\times 3N$ matrix $\widetilde{H}$, whose $(i,j)^{th}$ block of size $3\times 3$ is given by

[TABLE]

Then, $\widetilde{H}$ admits the decomposition in (8.5), as required. Thus, it only remains to show how to extract the set $\{\widetilde{s}_{ij}\}_{i<j\in[N]}$ from the values $\{\widetilde{s}_{ij}\widetilde{s}_{jk}\}_{i\neq j\neq k\in[N]}$ in (8.13).

Let us define the ${\binom{N}{2}}\times{\binom{N}{2}}$ matrix

[TABLE]

where $i<j\in[N]$, $k<l\in[N]$, and the products $\{\widetilde{s}_{ij}\widetilde{s}_{jk}\}_{i\neq j\neq k\in[N]}$ are computed using (8.13). The following proposition, the proof of which is given in Appendix A.8, shows that the signs $\widetilde{s}_{ij}$ can be extracted from $S$ in (8.16).

Proposition 8.2**.**

The leading eigenvalue of $S$ in (8.16) is $2(N-2)$ and it is simple. Moreover, define $u_{s}=(\widetilde{s}_{ij})_{i<j\in[N]}$ to be the vector of length $\binom{N}{2}$ with entries $\widetilde{s}_{ij}$. Then, $u_{s}$ is an eigenvector of $S$ corresponding to the eigenvalue $2(N-2)$.

By Proposition 8.2, the eigenvector $u_{s}$ of $S$ in (8.16) gives the set $\{\widetilde{s}_{ij}\}_{i<j\in[N]}$. The procedure for the signs adjustment (and the construction of $\widetilde{H}$ of (8.5)), is summarized in Algorithm 5.

9 Numerical experiments

We implemented Algorithms 1–5 in Matlab and tested them on a dataset of raw projection images of the beta-galactosidase enzyme [4], which has a $D_{2}$ symmetry. All tests were executed on a dual Intel Xeon E5-2683 CPU (32 cores in total), with 768GB of RAM running Linux, and four nVidia GTX TITAN XP GPU’s. Section 9.1 provides some of the implementation details for Algorithms 1–5, and Section 9.2 presents the results on the experimental dataset.

9.1 Implementation details

To execute Algorithm 1, we need to discretize the space of rotations $SO(3)$. To that end, we generated a pseudo-uniform spherical grid of $K$=1200 points $z_{k}$ on $S^{2}$, using the Saaf-Kuijlaars algorithm [10]. Then, for each $z_{k}=(a_{k},b_{k},c_{k})^{T}\in S^{2}$ on the spherical grid, we computed the set of rotations

[TABLE]

where $\theta_{l}=2\pi l/L$ for $l=0,1,\ldots,L-1$, and the vectors $u_{k}$ and $w_{k}$ are given by

[TABLE]

It is easily verified that the vectors $w_{k}$, $u_{k}$ and $z_{k}$ form an orthonormal set, and that $Q_{kl}\in SO(3)$. The third column $z_{k}$ of each rotation $Q_{kl}$ is the beaming direction corresponding to the rotation $Q_{kl}$, and the vectors

[TABLE]

are the coordinate systems for the plane perpendicular to $z_{k}$. Thus, each set of matrices $\{Q_{kl}\}_{l=0}^{L-1}$ where $k\in\{1,\ldots,K\}$ is a discretization of the set of rotations in $SO(3)$ with beaming direction $z_{k}$. We found experimentally that choosing $L=72$ (together with $K=1200$) is sufficient to obtain accurate results.

As for runtime, for a set of 500 projection images, it took 1512 seconds to compute all sets of relative rotations $\{R_{i}^{T}g_{m}R_{j}\}_{m=1}^{4}$ (Algorithm 1), 720 second to synchronize handedness (Algorithm 2), 5784 seconds to compute the partition in (7.1) (Algorithms 3 and 4), and 1378 seconds to adjust the signs $s_{ij}$ of (7.1) (Algorithm 5).

9.2 Beta-galactosidase experimental results

We applied Algorithms 1–5 to the EMPIAR-$10061$ dataset [1] from the EMPIAR archive [6]. The dataset consists of 41,123 raw particles images, each of size $768\times 768$ pixels, with pixel size of $0.3185$ Å. To generate class averages from this dataset, we used the ASPIRE software package [2] as follows. First, all images were phase-flipped (in order to remove the phase-reversals in the CTF), down-sampled to size $129\times 129$ pixels (hence with pixel size of $1.9$ Å), and normalized so that the noise in each image has zero mean and unit variance. We then split the images into two independent sets, each consisting of 20,560 particle images, and all subsequent processing was applied to each set independently.

We next used the class-averaging procedure in ASPIRE [2] to generate 2000 class averages from each of the two sets of particle images (using the EM-based class averaging algorithm in ASPIRE). A sample of these class averages is shown in Fig. 3. The input to our algorithm was 500 out of the 2000 class averages (by selecting every 4th image).

Next, we used the algorithms described in this paper to estimate the rotation matrices that correspond to the 500 class averages, and reconstructed the three-dimensional density map using the class averages and their estimated rotation matrices. The resolution of the reconstructed volume, assessed by comparing the reconstructions from the two independent sets of class averages is 8.23 Å, using the Fourier shell correlation (FSC) $0.143$-criterion [14] (Fig. 4(a)). When comparing our reconstructions to a high resolution reconstruction of the molecule (EMD-$7770$ [3]), the resolution estimated using the 0.5-criterion of the FSC is 9.88 Å (Fig. 4(b)).

10 Summary and future work

In this paper, we presented a procedure for estimating the orientations corresponding to a given set of projection images of a $D_{2}$-symmetric molecule. We have shown that the set of relative rotations between all pairs of images admits a special graph structure, and demonstrated that this structure can be exploited to recover the rotations. We then demonstrated our method by reconstructing an ab-initio model from an experimental set of cryo-EM images.

An obvious future work is to extend the proposed method to $D_{n}$ for $n\geq 3$. Preliminary theoretical analysis suggests that this can be achieved by combining the method of the current paper with the algorithms derived in [9].

Acknowledgments

This research was supported by the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (grant agreement 723991 - CRYOMATH).

Appendix A Appendix

A.1 Proof of Proposition 5.1

Let $(R_{i}^{T}g_{\tau(m)}R_{j})_{m=1}^{4}$ be a permutation of the 4-tuple $(R_{i}^{T}g_{m}R_{j})_{m=1}^{4}$.

This follows immediately from (5.2) and (5.3), by noting that the tuple $(\tau(2),\tau(3),\tau(4))$ is a permutation of $(2,3,4)$ . 2. 2.

Suppose that $\tau(2)=1$. Then, $(\tau(1),\tau(3),\tau(4))$ is a permutation of $(2,3,4)$. By (5.5) and (5.3), we have

[TABLE]

By (5.2) and (5.3), we have

[TABLE]

as required. The proof for the cases $m=3,4$ is similar.

∎

A.2 Proof of Proposition 6.1

First, we note that the set $\{g_{1},g_{2},g_{3},g_{4}\}$ of (1.3) forms a multiplicative group of matrices, known in literature as the Klein four-group. In particular, $g_{1}=I$ is the identity element of the group, and we have

[TABLE]

Now, since $R_{i},R_{j}$ and $R_{k}$ are in $SO(3)$, we have

[TABLE]

By (A.1), each member of the group $\{g_{m}\}_{m=1}^{4}$ is its own inverse, and so for any triplet $(m,l,r)\in\{1,2,3,4\}^{3}$ we have that $g_{\tau_{ij}(m)}g_{\tau_{jk}(l)}g_{\tau_{ki}(r)}=I$ if and only if

[TABLE]

Now observe, that since $\tau_{ij},\tau_{jk},\tau_{ki}\in S_{4}$ are permutations, we have

[TABLE]

from which it follows that there are 16 possible products $g_{\tau_{ij}(m)}g_{\tau_{jk}(l)}$ on the left-hand side of (A.3), corresponding to the 16 products $R_{ij}^{m}R_{jk}^{l}$ for $m,l\in\{1,2,3,4\}$. Thus, there are at most 16 triplets $(m,l,r)\in\{1,2,3,4\}^{3}$ for which (A.3) is satisfied. On the other hand, combining the closure property of groups with (A.4) gives us that there always exists an element $g_{\tau_{ki}(r)}$ of $\{g_{\tau_{ki}(m)}\}_{m=1}^{4}$ such that (A.3) is satisfied. We conclude that are exactly 16 triplets $(m,l,r)\in\{1,2,3,4\}^{3}$ such that (A.3) is satisfied, from which by (A.2), the proof is concluded.

∎

A.3 Proof of Theorem 7.2

We assume without loss of generality, that the first ${\binom{N}{2}}$ rows of $\Omega$ correspond to the vertices $\{v_{ij}^{1}\}_{i<j\in[N]}$, the following ${\binom{N}{2}}$ rows correspond to the vertices $\{v_{ij}^{2}\}_{i<j\in[N]}$, and the last ${\binom{N}{2}}$ rows correspond to $\{v_{ij}^{3}\}_{i<j\in[N]}$. Let $\Sigma^{+}$ be the ${\binom{N}{2}}\times{\binom{N}{2}}$ matrix given by

[TABLE]

In [8], it was shown that the spectrum of $\Sigma^{+}$ is given by

[TABLE]

We now show how to relate the spectrum of $\Omega$ to that of $\Sigma^{+}$.

Let $B_{rp}$ for $r,p=1,2,3$, be an ${\binom{N}{2}}\times{\binom{N}{2}}$ matrix which consists of rows $((r-1){\binom{N}{2}}+1),\ldots,r{\binom{N}{2}}$ and columns $((p-1){\binom{N}{2}}+1),\ldots,p{\binom{N}{2}}$ of the matrix $\Omega$. That is, we partition $\Omega$ into 9 block matrices of equal dimensions. By assumption, for any $r\in\{1,2,3\}$, both the rows and columns of the matrix $B_{rr}$ correspond to vertices of the same set $C_{r}$ in (7.3). Thus, by (7.4) we have that for $r\in\{1,2,3\}$,

[TABLE]

which by (A.5) gives that $B_{rr}=\Sigma^{+}$ for all $r\in\{1,2,3\}$. Now, consider any matrix $B_{rp}$ for $r\neq p\in\{1,2,3\}$, and note that its rows correspond to the vertices of the set $C_{r}=\{v_{ij}^{r}\}_{i<j\in[N]}$, whereas its columns correspond the vertices of the set $C_{p}=\{v_{ij}^{p}\}_{i<j\in[N]}$. Again, by (7.4), we have that for $r\neq p\in\{1,2,3\}$,

[TABLE]

Thus, by (A.5), we have that $B_{rp}=-\Sigma^{+}$ for all $r\neq p\in\{1,2,3\}$, and we conclude that

[TABLE]

Now, let $u\in\mathbb{R}^{\binom{N}{2}}$ be any eigenvector of $\Sigma^{+}$ corresponding to an eigenvalue $\lambda$. Denote $z=(0,\ldots,0)^{T}\in\mathbb{R}^{\binom{N}{2}}$, and consider the column vectors of length $3{\binom{N}{2}}$

[TABLE]

By (A.9), we have

[TABLE]

This shows that if $\lambda$ is an eigenvalue of $\Sigma^{+}$, then $-\lambda$ and $2\lambda$ are eigenvalues of $\Omega$. Furthermore, if $u$ is an eigenvector of $\Sigma^{+}$ with eigenvalue $\lambda$, then $u_{1}$ in (A.10) is an eigenvector of $\Omega$ corresponding to the eigenvalue $-\lambda$ of $\Omega$, and $u_{2}$ and $u_{3}$ in (A.10) are eigenvectors of $\Omega$ corresponding to the eigenvalue $2\lambda$ of $\Omega$. Now, note that

[TABLE]

and thus $u_{2}$ and $u_{3}$ in (A.10) are independent for any $u\in\mathbb{R}^{\binom{N}{2}}$. Furthermore, since all eigenvalues of $\Sigma^{+}$ are non-zero (see (A.6)), we have that $2\lambda\neq-\lambda\neq 0$, and so $u_{1}$ is in a different eigenspace of $\Omega$ than $u_{2}$ and $u_{3}$. Thus, all three vectors in (A.10) are independent eigenvectors of $\Omega$. Let us denote the multiplicity of an eigenvalue $\lambda$ of the matrix $\Sigma^{+}$ by $m_{\Sigma^{+}}(\lambda)$, and the multiplicity of an eigenvalue $\mu$ of $\Omega$, by $m_{\Omega}(\mu)$. We also denote the three eigenvalues of $\Sigma^{+}$ by $\lambda_{1},\lambda_{2}$ and $\lambda_{3}$. It is simple to verify, that if a pair of vectors $u$ and $v$ are independent, then so are the pairs of vectors $\{(u,u,u)^{T},(v,v,v)^{T}\}$, $\{(u,z,-u)^{T},(v,z,-v)^{T}\}$, and $\{(u,-2u,u)^{T},(v,-2v,v)^{T}\}$. Thus, by (A.11), if $\lambda$ is an eigenvalue of $\Sigma^{+}$, then the eigenvalues $2\lambda$ and $-\lambda$ of $\Omega$ satisfy

[TABLE]

which gives us that

[TABLE]

On the other hand, since $\Omega$ has dimensions $3{\binom{N}{2}}\times 3{\binom{N}{2}}$, we have

[TABLE]

by which we have that

[TABLE]

for otherwise, by (A.12) we would have a strong inequality in (A.13), which is a contradiction to (A.14).

We conclude that the set of eigenvalues of $\Omega$ is given by $\{2\lambda_{i},-\lambda_{i}\}_{i=1}^{3}$. Finally, the multiplicities of the eigenvalues of $\Omega$ in (7.5) are computed by combining (A.6) with (A.15).

∎

A.4 Proof of Proposition 7.4

We begin by introducing some notation and definitions which we use in the current and subsequent proofs. Let $P_{\sigma}$ denote the $3\times 3$ permutation matrix of $\sigma\in S_{3}$, i.e., $P_{\sigma}$ satisfies $P_{\sigma}v=(v(\sigma(1)),v(\sigma(2)),v(\sigma(3))^{T}$ for any $v\in\mathbb{R}^{3}$. Thus, using the notation introduced in (7.6) we can write

[TABLE]

where $u_{\alpha}$ and $u_{\beta}$ are given in Definition 7.3, and $u_{3c}$ and $u_{2c}$ are defined in (7.7).

Now, using the notation of (7.9), consider the block $\Omega_{(i,j)(j,k)}$ of $\Omega$ for some $i<j<k\in[N]$. If it were the case that $\sigma_{ij}$ and $\sigma_{jk}$ are the identity permutations, then by (7.9) and (7.4), we would have that

[TABLE]

However, in general, $\Omega_{(i,j)(j,k)}$ is the $3\times 3$ block of $\Omega$ which we get by taking the entries of $\Omega$ in the rows corresponding to the triplet of vertices $(v_{ij}^{\sigma_{ij}(1)},v_{ij}^{\sigma_{ij}(2)},v_{ij}^{\sigma_{ij}(3)})$, and columns corresponding to the triplet of vertices $(v_{jk}^{\sigma_{jk}(1)},v_{jk}^{\sigma_{jk}(2)},v_{jk}^{\sigma_{kl}(3)})$. In other words, $\Omega_{(i,j)(j,k)}$ is obtained from the matrix in (A.17) by permuting its rows by $\sigma_{ij}$ and its columns by $\sigma_{jk}$, that is

[TABLE]

By the same argument applied to $\Omega_{(i,j)(k,l)}$ whenever $|\{i,j\}\cap\{k,l\}|=1$, we have that

[TABLE]

The reason that $P_{\sigma_{kl}}$ in (A.18) is transposed, is that in order to permute the columns of a matrix by $\sigma_{kl}$, one has to multiply it on the right by $P_{\sigma_{kl}}^{T}$. We now prove Proposition 7.4

Proof of Proposition 7.4.

Let us compute the Rayleigh quotient for $u_{\alpha}$. By the second equality in (7.14), we have

[TABLE]

By (A.16) and (A.18), and since permutation matrices are orthogonal, we have

[TABLE]

It is straightforward to check that each term in the sum (A.20) equals exactly $24\alpha^{2}$. By Definition 7.3, $\alpha=(2{\binom{N}{2}})^{-\frac{1}{2}}$, and since there are $\binom{N}{3}$ triplets $i<j<k\in[N]$, the sum in (A.20) amounts to

[TABLE]

By Theorem 7.2, $\mu_{c}=4(N-2)$ is the leading eigenvalue of $\Omega$, and thus, $u_{\alpha}$ maximizes the Rayleigh quotient of $\Omega$, by which we have that $u_{\alpha}$ is in the eigenspace of $\mu_{c}$. A similar calculation for $u_{\beta}$ shows that it is also in the same eigenspace. Finally, observe that since $P_{\sigma_{ij}}$ are orthogonal for all $i<j\in[N]$, we have

[TABLE]

∎

A.5 Proof of Lemma 7.7

We begin by showing the first equality in (7.14). Fix some $i<j\in[N]$, and observe that for any $k<l\in[N]$ such that $|\{i,j\}\cap\{k,l\}|=1$, we have that

[TABLE]

where in the first case of (A.21) it cannot be that $i=k$ since then we would have that $\left\lvert\{i,j\}\cap\{k,l\}\right\rvert>1$ (and similarly for the other cases). From (A.21) and (7.12) we get that

[TABLE]

This shows that for fixed $i<j$ it holds that

[TABLE]

and so by taking a union over all $i<j\in[N]$, we get from (7.11) that

[TABLE]

Conversely, suppose that $(i,j)(k,j)\in A_{ij}^{1}$ for some $i<j\in[N]$ and $k<j$, $k\neq i$. Then, we have that $(i,j)(k,j)=(i,j)(k,l)$ for $l=j$, $|\{i,j\}\cap\{k,l\}|=1$ and $k<l$. Thus, by (7.11), for any $k\in[N]$ such that $k<j$ and $k\neq i$ we have that $(i,j)(k,j)\in A$, from which we have that $A_{ij}^{1}\subseteq A$. Applying a similar argument to $A_{ij}^{2},A_{ij}^{3}$ and $A_{ij}^{4}$, we get that $A_{ij}^{2},A_{ij}^{3},A_{ij}^{4}\subseteq A$ for all $i<j\in[N]$, from which we get that

[TABLE]

Taking a union over all $i<j\in[N]$ we have

[TABLE]

which together with (A.22) proves the first equality in (7.14).

Let us now show the second equality in (7.14). Suppose that $i<j<k\in[N]$. Then, from (7.12) we have that

[TABLE]

By (A.24) and the definition of $A_{ijk}$ and $A_{ijk}^{f}$ in (7.13), we have

[TABLE]

Thus, taking the union over all $i<j<k\in[N]$ gives us that

[TABLE]

Now, for each $i<j<k\in[N]$ we have

[TABLE]

by which we have that (renaming the indices to those given in (7.12))

[TABLE]

Taking a union over all $i<j<k\in[N]$ we get

[TABLE]

thus, by (A.25), we have

[TABLE]

where the last inequality follows from (A.23).

Conversely, take $(i,j)(k,l)\in A$. By the first part of the proof, $(i,j)(k,l)$ belongs to one of the sets in (7.12). If $(i,j)(k,j)\in A_{ij}^{1}$, then we have that either $k<i<j$ or $i<k<j$, and thus

[TABLE]

that is, $(i,j)(k,j)\in A_{kij}$ or $(i,j)(k,j)\in A_{ikj}^{f}$. This shows that (renaming the indices to the order given in (7.13))

[TABLE]

for either $k<i<j\in[N]$ or $i<k<j\in[N]$, by which we conclude that

[TABLE]

In the same manner one can show that

[TABLE]

for all $i<j\in[N]$. Thus, we have that

[TABLE]

for all $i<j\in[N]$. Taking a union over all $i<j\in[N]$ and using the first equality in (7.14) gives

[TABLE]

which concludes the proof of the second equality in (7.14).

∎

A.6 Proof of Proposition 7.8

Define the sets

[TABLE]

where $u_{2c}$ and $u_{3c}$ were defined in (7.7), and $P_{\sigma}$ is the permutation matrix of $\sigma\in S_{3}$. That is, $\Gamma_{3c}$ is the set of all permutations of $u_{3c}$ and $\Gamma_{2c}$ is the set of all permutations of $u_{2c}$.

Lemma A.1**.**

Let $u,v\in\mathbb{R}^{3}$ be such that $u=P_{\sigma}u_{3c}$ and $v=P_{\sigma}u_{2c}$ for some $P_{\sigma}\in S_{3}$. Suppose, that $u^{\prime},v^{\prime}\in\mathbb{R}^{3}$ are such that $(u^{\prime}\;v^{\prime})=(u\;v)R(\theta)$ for some $\theta\in[0,2\pi)$, $u^{\prime}\in\Gamma_{3c}$ and $v^{\prime}\in\Gamma_{2c}$, where $R(\theta)$ was defined in (7.28). Then, there exists $P_{\sigma}^{\prime}\in S_{3}$ such that

[TABLE]

Proof.

First, since $P_{\sigma}$ is orthogonal, then

[TABLE]

Define $W=\{(x,y,z)^{T}|\;x+y+z=0\}$, and note that $W$ is a linear subspace of $\mathbb{R}^{3}$ that contains $\Gamma_{3c}$ and $\Gamma_{2c}$. Since $(u^{\prime}\;v^{\prime})=(u\;v)R(\theta)$ and $R(\theta)$ is an orthogonal matrix, by (A.28) we have that $u^{\prime}$ and $v^{\prime}$ are also orthogonal vectors, and by assumption $u^{\prime}$ and $v^{\prime}$ are in $\Gamma_{3c}$ and $\Gamma_{2c}$, respectively, and so they are in $W$. Now, since $v,v^{\prime}\in\Gamma_{2c}$ there exists $P_{\sigma}^{\prime}\in S_{3}$ such that $v^{\prime}=P_{\sigma}^{\prime}v$. Since $P_{\sigma}^{\prime}$ is an orthogonal matrix, we have

[TABLE]

Finally, since $W$ is of dimension 2, there are exactly two vectors perpendicular to $v^{\prime}$ in $W$. Thus, it must be that $u^{\prime}=\pm P_{\sigma}^{\prime}u$, from which we get (A.27). ∎

We now prove Proposition 7.8.

Proof of Proposition 7.8.

The function $f_{c}(\theta)$, given in (7.35), satisfies $f_{c}(\theta)\geq 0$ for all $\theta\in[0,2\pi)$. Suppose that $\theta$ is such that $f_{c}(\theta)=0$. Such a $\theta$ necessarily exists, since if we choose $\theta$ in (7.33) such that $(v_{a}^{\theta},v_{b}^{\theta})=(\pm u_{\alpha},u_{\beta})$ (which can be done due to (7.32)), then (7.35) equals zero. In the notation of (7.6), we now show that there exists $\sigma\in S_{3}$ such that either

[TABLE]

or that

[TABLE]

from which it follows that $(v_{a}^{\theta},v_{b}^{\theta})\in\{\pm(u_{\alpha}^{\sigma},u_{\beta}^{\sigma})\;|\;\sigma\in S_{3}\}$.

First, we show that

[TABLE]

Indeed, for each pair $i<j\in[N]$, looking at the first square brackets in (7.35), we have that

[TABLE]

This shows that $(v_{a}^{\theta})_{ij}$ must be a permutation of $u_{3c}$ in (7.7), that is, $(v_{a}^{\theta})_{ij}\in\Gamma_{3c}$. Similarly, looking at the second square brackets in (7.35), we have that

[TABLE]

[TABLE]

which is possible only if $(v_{b}^{\theta})_{ij}$ is a permutation of the vector $u_{2c}$ in (7.7), i.e., $(v_{b}^{\theta})_{ij}\in\Gamma_{2c}$, which shows (A.31).

Now, by (7.32) and (7.33), the vectors $v_{a}^{\theta}$ and $v_{b}^{\theta}$ are either given by

[TABLE]

or by

[TABLE]

for some $\varphi\in[0,2\pi)$, where $v_{a}$ and $v_{b}$ is the pair of orthogonal eigenvectors of $\Omega$ defined in (7.25). Let us assume first the case (A.32). In the notation of (7.6), we define

[TABLE]

Also, denoting by

[TABLE]

the $2\times 3$ matrices with columns $\pm u_{3c}$ and $u_{2c}$ (defined in (7.7)), we get that by Definition 7.3 and (7.6), for each pair $i<j\in[N]$ we have

[TABLE]

In particular, for $i=1$ and $j=2$ we have

[TABLE]

Thus, by (A.35) and (A.36) it follows that for all $i<j\in[N]$ we have

[TABLE]

Now, by (A.32) and (7.6) we have that

[TABLE]

By (A.36) we have that $(u_{\alpha})_{12}=P_{\sigma_{12}}u_{3c}$ and $(u_{\beta})_{12}=P_{\sigma_{12}}u_{2c}$, and by (A.31), we have that $(v_{a}^{\theta})_{12}\in\Gamma_{3c}$ and $(v_{b}^{\theta})_{12}\in\Gamma_{2c}$. Thus, using (A.38), by Lemma A.1 there exists a permutation matrix $P_{\tau}$ such that either

[TABLE]

First, assume that the case on the left of (A.39) holds. It then follows from (7.6), (A.32), (A.37) and (A.36) that

[TABLE]

Writing $P_{\sigma}=P_{\sigma_{12}}^{T}P_{\tau}P_{\sigma_{12}}$, we get by (A.34) that for all $i<j\in[N]$

[TABLE]

By (7.8) and (7.6), we have that

[TABLE]

for all $i<j\in[N]$. The last two equations show that (A.29) holds, which proves the proposition for the case $\eqref{eq:Prop8Case1}$, when the identity on the left of (A.39) holds.

If $\eqref{eq:Prop8Case1}$ holds as well as the case on the right of (A.39), where $A_{\alpha\beta}^{+}R(\theta+\varphi)=P_{\tau}A_{\alpha\beta}^{-}$, then by repeating the latter calculation with $A_{\alpha\beta}^{+}$ and $A_{c}^{+}$ replaced by $A_{\alpha\beta}^{-}$ and $A_{c}^{-}$, we get that

[TABLE]

for all $i<j\in[N]$, i.e., that (A.30) holds, which proves the proposition for the case $\eqref{eq:Prop8Case1}$ when the identity on the right of (A.39) holds. This concludes the proof for the case $\eqref{eq:Prop8Case1}$.

In the case where $\eqref{eq:Prop8Case2}$ holds, we get by the same method of proof that either (A.30) or (A.29) hold, which proves the proposition for this case.

To conclude, we have shown that given $\theta$ which minimizes (7.35), the vectors $v_{a}^{\theta}$ and $v_{b}^{\theta}$ defined in (7.33), must satisfy either $(v_{a}^{\theta},v_{b}^{\theta})=(u_{\alpha}^{\sigma},u_{\beta}^{\sigma})$ or $(v_{a}^{\theta},v_{b}^{\theta})=(-u_{\alpha}^{\sigma},u_{\beta}^{\sigma})$ for some $\sigma\in S_{3}$.

∎

A.7 Proof of Proposition 8.1

Suppose that (8.6) holds, and fix an arbitrary $n\in[N]$. Then, the $(i,j)^{th}$ $3\times 3$ block of $H$ is given by $s_{ij}v_{i}^{T}v_{j}=s_{in}s_{nj}v_{i}^{T}v_{j}$. Thus, since $s_{nj}=s_{jn}$ for all $j\in[N]$, we have

[TABLE]

which gives (8.7) and shows that $H$ indeed has rank 1.

Now suppose that (8.6) does not hold, and assume without loss of generality that $s_{12}s_{23}=-s_{13}$. Denote the $m^{th}$ entry of a row $v_{i}$ by $v_{i}(m)$, $i\in[N]$. The rank 1 matrix $v_{1}^{T}v_{2}$ is non-zero, thus, there exist $r,l\in\{1,2,3\}$ such that $v_{1}(r),v_{2}(l)\neq 0$, that is, the $r^{th}$ and $l^{th}$ entries of the vectors $v_{1}$ and $v_{2}$, respectively, are non-zero. By (8.2), we have that $s_{11}=1$, and thus, the first three rows of $H$ are given by the $3\times 3N$ matrix

[TABLE]

Similarly, the next three rows of $H$ are given by the $3\times 3N$ matrix

[TABLE]

Thus, since each vector $v_{i}$ is of length 3, rows number $r$ and $3+l$ of $H$ are given by

[TABLE]

Multiplying (A.44) by $\frac{1}{v_{1}(r)}$ and (A.45) by $\frac{s_{12}}{v_{2}(l)}$, by our assumption we get

[TABLE]

Since $v_{3}\neq 0$, the latter two vectors are linearly dependent only if $s_{13}=-s_{13}$, which is impossible. Therefore, rows number $r$ and $3+l$ of $H$ given in (A.44) and (A.45) are linearly independent, which implies that $\operatorname{rank}(H)\geq 2$.

∎

A.8 Proof of Proposition 8.2

Fix a pair of indices $i<j\in[N]$. We begin by deriving an expression for $(Su_{s})_{ij}$, which is the entry of $Su_{s}$ corresponding to the $(i,j)^{th}$ row of $S$. By (8.16) and the first equality in (7.14) of Lemma 7.7, we have

[TABLE]

where the 4 sums to the right of the first equality in (A.46) correspond to the 4 sets $A_{ij}^{m}$ of Lemma 7.7. Since there are exactly $N-2$ indices $k\in[N]$ for which $k\neq i,j$, we conclude from (A.46) that $(Su_{s})_{(i,j)}=2(N-2)\widetilde{s}_{ij}$ for all $i<j\in[N]$. Thus, $u_{s}$ is an eigenvector of $S$ corresponding to the eigenvalue $2(N-2)$. Next, we show that $2(N-2)$ is simple. Define a partition of the set $\{\widetilde{s}_{ij}\}_{i<j\in[N]}$ of (8.12) into two disjoint sets

[TABLE]

and note that a pair $\widetilde{s}_{ij}$ and $\widetilde{s}_{kl}$ are in the same set of (A.47) iff $\widetilde{s}_{ij}\widetilde{s}_{kl}=1$. Thus, by (8.16) and (A.47), the matrix $S$ is given by

[TABLE]

In [8], it was shown that the leading eigenvalue of $S$ is simple and is given by $2(N-2)$, which concludes the proof.

∎