Double-tilt electron tomography

ultramicroscopy ELSEVIER

Ultramicroscopy 60 (1995) 393-410

Double-tilt electron tomography Pawel Penczek, Michael Marko, Karolyn Buttle, Joachim Frank

*

Wadsworth Center for Laboratories and Research. New York State Department of Health, P.O. Box 509, Empire State Plaza, Albany, NY 12201-0509. USA

Received 14 March 1995; accepted 10 June 1995

Abstract Fidelity of tomographic reconstructions is improved and reconstruction artifacts are reduced, without increasing the number of projections, by combining tilt series taken around two orthogonal axes. Test reconstructions were made from high-voltage EM of rat liver mitochondria in a 0.6 pm thick plastic section. A number of schemes for selecting tilt angles for the projections are compared. A new method for aligning fiducial markers is described. It uses an iterative algorithm to determine the shift, scale, in-plane rotation and tilt angle for each tilt image, enforcing agreement of the expected locations of the fiducial markers in 3D space. These 3D locations are used to find the orientation between two tilt series and to merge both sets of projections.

1. Introduction

structural biology laboratories, and its use is certain to expand, especially when applied to the thicker specimens that can be studied with high- and intermediate-voltage electron microscopes, possibly combined with energy filtering. For objects in the size range of many cellular organelles, no other method can provide better 3D resolution. New findings from familiar material such as mitochondria 111, kinetochores [2], chromatin and chromosomes [3-51, synaptonemal complex [6] and Golgi apparatus [7], challenge standard textbook conceptions. In the quest for higher fidelity of electron tomographic reconstructions, numerous technical as well as conceptual limitations have to be addressed and, if possible, overcome. These include radiation damage, section shrinkage and angular restrictions. The latter arise, when sections are used, from the fact that the section has a surface area much larger than its thickness; thus in one direction it becomes practically impenetrable for the electron beam. This, in turn, results in artifacts in the reconstructed volume, which are caused by missing high-tilt projection data. In the current work, we put forward a method that overcomes the problem of missing data at least partially. Electron

l

tomography

has become

a routine

tool in several

Corresponding author.

0304-3991/95/$09.50

0 1995 Elsevier Science B.V. All rights reserved

SSDI 0304-3991(95)00078-X

394

P. Penczek et al./ Ultramicroscopy 60 (1995) 393-410

A second single-tilt series is collected after rotation of the specimen by ninety degrees around an axis perpendicular to its support plane [8]. This results in an increase of Fourier information collected and, owing to the resulting more-even distribution, in a reduction of reconstruction artifacts. Two computer programs were developed, one designed for the purpose of alignment of tilt images using fiducial markers and another for subsequent determination of rotation angles between the two tilt series. These programs are both described in the Appendix. The addition of a second single-tilt series implies not only increased time and effort to collect more data, but, more importantly, increased radiation damage to the specimen. In an effort to overcome these drawbacks, and yet optimize the quality of the reconstruction, we investigated a number of tilt-image selection schemes, in which we used a reduced set of projections. Taking into account the particular geometry of double-tilt reconstruction, it would be expected that high-tilt images should contribute relatively more to the quality of the 3D reconstruction. (Excessive specimen thickness for the acceleration voltage used can diminish this gain due to the increased proportion of inelastically scattered electrons and the resulting blur.) This hypothesis was tested in two ways. The tilt angle increment was adjusted so that more projections with higher tilt than projections with lower tilt were used. Alternatively, after one complete single-tilt series was already collected, only a number of high-tilt projections from the second series was added. In the former case we found that it is possible to obtain a more faithful reconstruction without using a greater number of tilt images than are used in a conventional single-tilt reconstruction.

2. Methods

2.1. Background 2.1 .I. Factors influencing reconstruction

quality The reconstruction quality of biological specimens using electron tomography cannot be properly discussed without considering the problem of missing Fourier information. Because of tilt-stage and specimen-thickness limitations, sectioned specimens are rarely tilted beyond 70”, and a maximum tilt of 60” is common with most standard goniometers. This means that a large part of the Fourier information remains unknown, which results in various distortions of the reconstructed object. This is known as the missing wedge problem. A similar problem in the random conical tilt scheme of data collection used in single-particle reconstruction [9] is referred to as the missing cone problem. Both terms refer to the shape of the missing regions in Fourier space. In both cases, much attention has been paid to the nature of the distortions due to the missing information, and their effects on the reconstruction [9-l I]. These distortions can be characterized (using a first-order approximation) by the elongation of the main lobe of the point-spread function. This causes a “blur” in the direction of missing region in the reconstructed object. Under some circumstances, a more important problem is the “weakening” (or virtual disappearance) of features that should have corresponding Fourier information in the missing region. This effect can be potentially more damaging than the blur, since it is strongly object-dependent and thus more difficult to analyze and predict. Comparison of conical-tilt and single-tilt data collection geometries [9] leads to the conclusion that, for the same maximum tilt angle, the conical geometry provides more Fourier information. However, the conical tilt method is difficult to carry out using commonly available tilt-rotation stages because the object cannot be maintained at the eucentric point. A substantial improvement in reconstruction quality should still be possible to achieve by collecting the projeetion data as two orthogonal single-axis tilt series. The change of axis can be accomplished by rotating the specimen by 90”, either by making use of the rotation capability of a tilt-rotation stage, or outside of the microscope when only a single-tilt stage is available. A double-tilt reconstruction will have the missing wedge in Fourier space reduced to a missing pyramid, which for high maximum tilt angle favorably compares with a missing cone at the same tilt angle.

P. Penczek et al./ Ultramicroscopy 60 (19951393-410

395

2.1.2. Alignment using jiducial markers Fiducial markers, usually gold particles from 5 to 40 nm in diameter, are placed on one or both surfaces of a sectioned specimen. They provide a convenient set of coordinates for the numerical alignment procedure. All the necessary calculations are done using these coordinates, avoiding the need to do any image processing for alignment purposes. In addition, and of great importance as a check, the results provide accurate and detailed information about errors and irreducible distortions present in the input data. This information is, to a certain extent, independent of the object itself and provides an independent assessment of the resolution that can be expected from a given data set. The calculated coordinates of markers in 3D space are also used to find the orientation between two sets of projections of the same object in the double-tilt collection scheme (see Appendix). 2.2. Materials and methods 2.2.1. Specimen and microscopy The test specimen was a 0.6 pm thick plastic section of isolated rat liver mitochondria [l]. The section was collected on a folding grid without supporting film, and was stained by immersion in 2% aqueous uranyl acetate (60°C for 1 h) and Reynold’s lead citrate (20°C for 20 min). 15 nm colloidal gold was applied to one surface of the section. The images were recorded at a magnification of 25 000 X on DuPont Lo-Dose Mammography film [ 121 using the Albany AEI EM7 MkII high-voltage electron microscope operating at 1000 kV. Sixty-eight images were recorded at 2” intervals, from - 70” to + 64” around the y-axis, and another sixty-seven images at 2” intervals, from - 70” to + 62” around the x-axis. To accomplish the switch of tilt axis, the specimen rod was removed from the microscope column, and, after rotating the grid by 90”, reinserted. The goniometer used was non-eucentric. During the data collection, the dose rate was kept constant at about 10 e- A-’ s- ‘. Before recording the images, initial radiation damage was allowed to run to saturation (20 min), so that no more shrinkage during data collection would occur [13]. To limit further distortion of the specimen during focusing and tilting, a low-light-level TV camera and an image wobbler were used. The exposure time was proportional to the absolute value of the tilt angle. Each series required about 1.5 h of beam time to complete. 2.2.2. Image digitization and marker picking Negatives were digitized using a modified Dage Model 81 camera and Androx ICS-400 image-processing board [ 141 in a Sun workstation. The images were approximately aligned in real time during digitization, making use of selected gold particles, with the aid of an axis-centered film stage and a camera monitor. An image size of 940 X 940 pixels and a pixel size of 3 nm were used. All the images were normalized such that the mass of the specimen remained constant for all the tilt angles. The coordinates of selected gold particles were found on each image by manual use of a cursor on the monitor (for this work we used twelve particles). These markers were selected such that they were evenly distributed around the object, so that the centers of gravity of both the set of markers and the object itself roughly coincided. This facilitated further manipulation of the data set; in particular, the determination of the rotation of the second tilt series. Each particle was given the same number in both series. The position of these particles was refined within each series by using a peak-search algorithm within a small box centered on each manually picked location. 2.2.3. Alignment and 30 reconstruction The numerical alignment method for an individual single-tilt series is described in Appendix A.1. Its application results in two sets of aligned tilt-series images and a pair of 3D coordinates (one from each tilt series) for each of the chosen gold particles on the specimen.

396

P. Penczek et al./ Ultramicroscopy 60 (1995) 393-410

Once the rotation of the second set of markers is established, the three Eulerian angles are used to calculate the angles defining the orientations of the micrographs in the second tilt series. The images themselves do not need to be altered. Then, a “merged” 3D reconstruction is calculated, using either the general-weighted back-projection algorithm [15] or iterative reconstruction algorithms [16]. Because of the large number of 3D reconstructions carried out in this study, the former method was used on account of its speed. For the same reason, we compressed the projection data. The aligned images were windowed for the smallest field occupied by the object (640 X 640 pixels), and then reduced in size to 320 X 320 pixels, each pixel representing the average of four neighboring pixels in the original image. This resulted in a 6 nm voxel size in the reconstructed volumes. The full-size 3D reconstruction with 3 nm voxel size was also calculated and will be presented in a forthcoming paper [Mannella et al., manuscript in preparation]. To assess the quality of the double-tilt 3D reconstructions we calculated, additionally, two single-axis tilt 3D reconstructions, the first one from the y-tilt set, the second one from the x-tilt set. Both single-axis tilt reconstructions were calculated in such a way that the orientation of the resulting volumes coincided with the orientation of the double-tilt reconstruction. 2.2.4. Tilt image selection In addition to comparisons with single-axis tilt reconstructions, some tests were performed to investigate whether the number of projections required for the double-tilt reconstruction can be lowered without significant loss of the reconstruction quality. Four test reconstructions were calculated using different subsets of the total set of 135 projection images available. The first test was a double-tilt reconstruction calculated using only tilt images below 60”. The purpose was to evaluate whether images with tilt angles higher than this value contribute significantly to the reduction of artifacts. Since many microscope stages lack the ability to tilt above 60” this question is of considerable practical importance. The second test reconstruction was calculated using evenly distributed projections in 4” steps (using every second image from the original 135image data set), reducing the total number of images to 68. The third and fourth test reconstructions were calculated to explore uneven projection-data distribution schemes (alternative projection distribution schemes for thin sections in single-tilt axis reconstructions were tested by Saxton et al. 1171 and, independently, by Levy et al. [18]). In the third reconstruction, to test the assumption that high-tilt images are relatively more important for the reduction of artifacts, we used a graduated tilt-angle increment for both single-axis tilt series in the following way: an 8” interval in the tilt range - 28” to + 28”; a 4” interval in the tilt ranges - 36” to - 48” and + 36” to f48”; and a 2” interval in the highest tilt ranges - 52” to - 70” and + 52” to + 64” (+62” for the second series). The total number of images was 65. The fourth reconstruction was made to test whether an already computed single-axis reconstruction can be improved by recording a small number of additional high-tilt images after rotating the grid by 90”, assuming that the specimen is still found to be in good condition. Thus, we used the full 70” tilt range around the y-axis, with the addition of tilt images between 50” and 70” around the x-axis. The total number of images was 86.

3. Results 3.1. Marker error analysis The set of twelve markers selected were aligned separately for each single-tilt series using a marker alignment procedure described in Appendix A.l. The results are shown in Table 1 for the first and second series, respectively. The average mean-squares error in the first series was 1.5 pixels and that for the second series was 1.7 pixels. (Note that the error in pixels refers to the original resolution of the 940 X 940 pixel

397

P. Penczek et al./ Ultramicroscopy 60 (1995) 393-410

Table 1 Results of alignment error of each point) Point

of markers (twelve points were aligned and their estimated coordinates

Mean-square error

Estimated coordinates Z

Y

X

(al Results of thefirst single-axis tilt series containing 68projections 313.9 15.7 1 189.8 - 209.1 2 11.6 - 404.0 3 - 149.5 - 346.4 4 - 269.2 - 269.1 5 -313.8 - 102.1 6 - 290.2 - 39.8 7 -312.2 260.3 8 43.0 382.6 9 165.8 281.1 10 296.6 231.7 11 314.2 200.3 12 Average error

in 3D space are given as well as mean-square

48.6 -7.1 - 64.4 - 79.7 - 88.7 -67.7 - 54.3 -9.2 77.6 72.7 88.3 84.0

1.2 1.4 1.7 1.2 1.1 1.4 1.3 1.8 1.5 2.0 1.4 1.3 1.5

1.8 1.7 1.9 1.3 1.3 1.8 1.8 2.7 1.5 1.3 1.5 1.4 1.7

Cb)Results of the second single-axis tilt series containing 64 projections 1 2 3 4 5 6 I 8 9 10 11 12 Average error

-

306.0 232.2 98.0 - 73.6

- 87.9 161.8 393.7 373.5

27.3 -8.1 - 47.4 - 52.8

208.3 287.6 277.7 363.3

325.8 171.1 104.0 - 183.8

- 56.7 - 39.4 - 29.3 0.4

- 38.8 102.7 243.6 266.8

- 384.1 -312.2 - 294.5 - 267.3

53.3 46.8 55.5 50.5

images.) The differences among errors of individual markers were not large. For the first series the errors were between 1.1 and 2.0; for the second between 1.3 and 2.7. No single marker had a particularly large error. A decreased average error could have been obtained by elimination of some markers, but this improvement was thought likely to be local, while the remaining areas of the object might suffer from increased distortions. It is worth noting that, since all the markers selected were located outside of the structure, the error within the field of interest is probably smaller than the alignment error for the markers. We assume that the observed error is caused by distortions of the specimen during irradiation, not by continued mass loss. The three-dimensional coordinates for both sets of markers were aligned to one another using a procedure described in Appendix A.2. The Eulerian angles found were + = - 46.9, 8 = 4.9 and cp= -55.8, with an additional scale adjustment for the second series of 1.002. The values of the angles indicate that the grid was tilted by 4.9” when it was repositioned in the stage and that the in-plane rotation of the grid exceeded 90”. The errors for alignment of the two series are included in Table 2. The average error was 2.2 pixels. This error is the combined effect of the individual tilt-series alignment errors and a slight electron-optical distortion (data not shown).

398

P. Penczek et al./ Ultramicroscopy

Table 2 Results of 3D alignment of markers (estimated coordinates Eulerian angles and scale adjustment found) Point

I

2 3 4 5 6 7 8 9 IO 11 12 Total error =

3.2. Three-dimensional

60 (1995) 393-410

of second set of markers after 3D alignment

Estimated coordinates

and application

of the rotation using

Error

X

Y

2

16.6 - 208.9 - 403.5 - 345.5 - 269.0 - 101.3 - 38.5 259.8 380.8 279.5 230.7 199.3

315.5 191.7 13.1 - 149.4 - 270. I -315.0 -291.3 -314.2 42.1 164.9 297.8 314.9

50.8 - 4.6 - 66.4 - 80.6 - 89.4 -67.1 - 52.2 - 9.4 75.1 72.4 88.3 83.0

2.9 3.2 2.5 1.3 1.3 1.5 2.8 2.0 3.3 1.8 1.6 1.6 2.2

reconstructions

The quality of the reconstructions was judged by visual appearance. In detail, we examined (1) the benefits of adding data from a second complete tilt series to a single-tilt set of projections, and (2) the viability of schemes that keep the number of projections the same as in single-tilt reconstruction. 3.2.1. Double-tilt reconstructions The complete double-tilt reconstruction using a maximum tilt angle of 70” is shown in Fig. 1. The y-axis and x-axis single-tilt reconstructions are shown in Fig. 2 and Fig. 3, respectively. Identical z-, y- and x-slices are

Fig. 1. Slices from a complete double-tilt reconstruction using a maximum tilt angle of 70”. The positions of the image of the x-slice (x) and the y-slice (y) within the volume are indicated by lines in the image of the z-slice (2). The diameter of the mitochondrion is 1.4 pm. Note good continuity of inner and outer mitochondrial membranes in slices along all three axes. In (y) and (x), the membranes are not seen at the bottom of the images because the specimen is cut by the surface of the section. In (xl, note the improved contrast and detail in the internal membrane, compared with Fig. 2x. Arrows: tubular connections between the intracrystal space and the intermembrane space.

P. Penczek et al. / Ultramicroscopy

60 (1995) 393-410

399

Fig. 2. Slices from a y-axis single-tilt reconstruction. Same slice positions as in Fig. 1. Note lack of membrane continuity in (y) (due the effects of the Fourier missing wedge) compared with Fig. ly.

shown for all volumes. The position of the X- and y-slices within the volume are indicated by lines in Fig. lz. The overall result is that features in the double-tilt reconstruction (Fig. 1) present a combination of the features seen in the single-tilt reconstructions separately. As expected, there is no significant improvement in the quality of the z-slices, since they are least affected by the distortions. The largest improvement can be observed on slices perpendicular to the tilt axes. For these slices (Fig. 2y and Fig. 3x), the effects of the Fourier “missing wedge” are most apparent, and a comparison with corresponding slices obtained for the double-tilt reconstruction (Fig. ly and Fig. lx, respectively) is most revealing. It can be seen that in Fig. 2y and Fig. 3x the membrane continuity at the upper edge of the slice is lost (note that the bottom edge does not show membranes because the section plane cuts the structure there). In the merged reconstruction (both Fig. ly and Fig. lx), the membrane continuity is preserved. In addition, the dark artifacts in Fig. 3y (a streaking effect due to the missing wedge in the reconstruction of nearby gold markers) have become strongly reduced in Fig. 1y. Interestingly,

Fig. 3. Slices from an x-axis single-tilt reconstruction. Same slice positions as in Fig. 1. Note lack of membrane continuity in (x) compared with Fig. lx. Jn (y) note dark streaking artifacts along upper edge of image (due the missing wedge in the reconstruction of nearby gold markers).

400

P. Penczek et al./ Ultramicroscopy

Fig. 4. Slices from a double-tilt nearly as good as Fig. 1.

reconstruction

using a maximum

60 (19951393-410

tilt angle of 60”. Same slice positions as in Fig. 1. Note that the quality is

some improvement occurs over the single-axis tilt reconstruction slices in the directions not affected by the missing wedge. Most notably, the internal membrane visible at the center of Fig. 2x appears to have better contrast and more detail in Fig. lx. A double-tilt reconstruction from the same data, but without tilt images above 60”, is shown in Fig. 4. Surprisingly, there is no significant difference between Fig. 4 and Fig. 1. The minor differences between two volumes are at the resolution limit and do not affect overall quality of the reconstruction. 3.2.2. Selection A double-tilt images (Fig. 5) A double-tilt

of a reduced

set of double-tilt images

reconstruction made with half of the total tilt images, simply by taking only the odd-numbered shows considerable deterioration of image quality, which is visible even in the z-slices. reconstruction made using the graduated tilt angle increment selection scheme described in

Fig. 5. Slices from a double-tilt reconstruction made with half of the total number of tilt images, by selectively taking only the odd-numbered images. Same slice positions as in Fig. 1. Note deterioration in quality along all three axes compared with Fig. 1.

P. Penczek er al. / Ullramicroscopy 60 (1995) 393-410

Fig. 6. Slices from a double-tilt reconstruction 1. Note improvement in quality (particularly membrane continuity.

401

made using a graduated tilt-angle increment selection scheme. Same slice positions as in Fig. in x- and y-slices) compared to Fig. 5, approaching the quality of Fig. 1 with respect to

2.2.4 is shown in Fig. 6. Visual inspection confirms the assumption that the higher-tilt images contribute more to the quality of the reconstruction that the lower-tilt images. Note that the membrane continuity in y- and x-slices is considerably better than in Fig. 5, which shows a reconstruction with the same number of projections, but with a constant tilt increment. In comparison with single-tilt reconstructions (Fig. 2 and Fig. 3), the structural information is more reliable even though the images do not appear as sharp, and membrane continuity is preserved in both y- and x-slices. Fig. 7 shows a reconstruction made from the full 70” tilt range around the y-axis, with the addition of tilt images between 50” and 70” around the x-axis. Contrary to expectations, this particular scheme failed to produce significant improvement over the corresponding single-axis tilt reconstruction (Fig. 2, in particular Section

Fig. 7. Slices from a double-tilt reconstruction made from the full 70” tilt range around the y-axis, with the addition of tilt images between 50” and 70” around the x-axis. Same slice positions as in Fig. 1. No significant improvement compared with Fig. 2.

402

P. Penczek

et al./

Ultramicroscopy

compare Fig. 2y and Fig. 7~1. As with the double-tilt relatively little contribution from the high-tilt data.

60 (19951393-410

reconstruction

without data above 60”, there appears to be

4. Discussion The analysis of marker alignment errors indicates that the behavior of the specimen in the microscope was fairly stable in the experiment performed. The average mean-square error within each of the series was found to be below two pixels, which indicates a high degree of alignment accuracy. For the purpose of a resolution estimation, we assume that the residual error has a Gaussian distribution since it results from the addition of many independent random errors. Equating resolution with the half-width of a Gaussian function, we can estimate that for a pixel size of 3 nm it should be possible to realize a resolution better than l/16 nm- ’ using the current method of recording HVEM data. Moreover, the relatively small average error between the two series (slightly exceeding two pixels) further confirms reproducibility of the results and proves that merging two sets of single-axis projection data is an effective method for reduction of artifacts caused by missing Fourier information. The improved quality of the merged reconstruction was illustrated in comparison with the two single-axis tilt reconstructions using projections with tilt angles extending to 70”. Imposing a limit of 60” to the tilt data in the merged reconstruction did not visibly deteriorate the reconstruction quality. Thus, we can suggest that double-tilt reconstructions, even without tilts beyond the 60” limit imposed by many microscope stages, should provide reconstructions superior to single-axis tilt reconstructions with the same tilt limits. Test reconstructions were performed to investigate alternative schemes of angular tilt data distribution in order to reduce the number of projections. The most obvious scheme, namely the use of an increased angular gap between projections in a double-tilt approach, was found not to be a good solution. Two schemes were designed to test the possibility of uneven angular distribution of projection data. On the basis of an analysis of all these schemes, we conclude that, in cases where the specimen can tolerate the additional irradiation, recording of a double-tilt series with a maximum tilt angle of 60” will suffice. Alternatively, in cases when the radiation dose is of concern (and thus the number of projections should be kept as small as possible) a graduated-angle double-tilt series should be collected. In the latter case a significant improvement of the reconstruction quality was achieved without the need to collect more projection data. This is of particular importance in future applications of the method proposed for tomography of ice-embedded specimens 119,201.

5. Conclusions The new approach to single-axis tilt 3D reconstruction presented in this paper helps to alleviate one of the major problems hindering further progress in high-resolution electron tomography: distortions caused by missing angular information. Double-tilt tomography visibly improves fidelity of the reconstruction by filling the Fourier space more evenly and reducing the “missing wedge” to a “missing pyramid”. The procedure proposed requires no modification of existing electron microscopes, and the additional effort involved in collection of the second single-axis tilt series remains within reasonable limits. The alignment algorithms (presented in the Appendix) are easy to implement and have been integrated (together with 3D tomographic software) into the SPIDER image processing system 1211. Preliminary examination of the test volumes reveals that a wealth of biological information can be expected to emerge from the double-tilt reconstruction. In our earlier work on the mitochondria [l], using single-axis tilt reconstruction, questions about the nature of the tubular connections (Fig. 1, arrow) between the intracristal space and the intermembrane space were raised. The much-improved visualization of the inner and outer membranes now allows study of membrane interactions in the vicinity of the entrance to the cristae.

P. Penczek et al./ Ultramicroscopy 60 (1995) 393-410

403

Acknowledgements We thank Dr. C.A. Mannella and Dr. B.F. McEwen for advice and suggestions. This work was supported, in part, by NIH PHS grant R01219 (Principal Investigator Dr. C. Rieder) and NSF grant BIR 921 9043 (Principal Investigator Dr. J. Frank), supporting the Biological Microscopy and Image Reconstruction Resource, NIH grant ROI GM 20169 (Principal Investigator Dr. J. Frank), and NSF grant MCB 921 9353 (Principal Investigator Dr. C.A. Mannella).

Appendix

A. Marker alignment

methods (P. Penczek and R. Feuerbach)

A.I. Marker alignment algorithm The fiducial-marker alignment algorithm, as originally developed by Lawrence [22] (see also Ref. [23]) and independently by Berriman et al. [24], seeks to minimize the least-squares error between the measured marker coordinates on the micrographs and the coordinates of the markers reprojected from a 3D marker model. The coordinates are calculated using transformations for rotation of the projection plane in 3D space (taking into account the three Eulerian angles), for in-plane translation of images, and for magnification changes. In its most general form, the resulting equation is non-linear. However, since the equation is differentiable, the minimization can be attempted by one of the available numerical algorithms, for example conjugate gradients. Even though in principle a minimum error can be found, the equation is not only non-linear, but the parameters to be determined are also highly correlated. This makes the task extremely difficult. Thus, motivated by the need for a fast and robust alignment algorithm using markers, we developed another approach. Using a number of simplifying assumptions (justified by the practice of single-axis tomography), we found that each of the required parameters, taken separately, can be estimated using closed-form solutions (i.e., not requiring iterations). Using these solutions, a marker-based alignment procedure is proposed. Marker coordinates in 3D space are estimated, and, using explicit solutions for each parameter, these coordinates are iteratively corrected. The algorithm is terminated when the projections of the 3D marker coordinates match the corrected coordinates of the fiducial markers on the micrographs with sufficient accuracy. The algorithm proves to be very fast and robust, and it does not require all the markers to be present on all micrographs. The estimated marker positions in 3D space for each data set separately were used to find their orientation using a least-squares, closed-form solution (originally suggested by Horn [251) (Appendix A.2). This enabled calculation of a “merged” 3D reconstruction. The purpose of the marker alignment algorithm is to determine a common coordinate system for a set of 2D projections, given the locations of fiducial markers in these projections. A set of parameters describing the orientation of projection in space includes Eulerian angles, in-plane translation and, to account for possible distortions, change of scale. In the following we assume that the microscope optical axis (beam direction) is parallel to the z-axis and the goniometer tilting axis coincides with the y-axis. For the purpose of simplifying the algorithm we will further assume that two Eulerian angles suffice to describe the orientation of a projection in space, which means that we assume that the goniometer tilting axis is strictly perpendicular to the microscope optical axis. Even though this assumption may not be strictly true in practice, any deviations are sufficiently small to justify it. Most importantly, this assumption reduces any dependencies between parameters and improves the stability of the algorithm. Furthermore, we assume that scale changes are isotropic; thus we have only one scale parameter per projection. It is also assumed that the initial tilt angles are known (they are read from the goniometer) and that initial orientations of projections are not far from the correct ones (which is true, provided that all the projections were scanned preserving original orientation). The set of parameters to be found has to correspond to the least-squares solution of the following set of equations: Pix’=

R,s,p,’ + di,

(1)

404

P. Penczek et al./ Ultramicroscopy

where i= 1, . . . . V, with v the number markers on each micrograph,

1

cos 0

Pi =

0

ei

is the projection

j = 1, . . . , s, with

(micrographs),

s the number

of

1

- sin

1

of projections

60 (1995) 393-410

0 61i

matrix for the ith tilt angle di,

Xj

[I

xj=

yj Zj

is the coordinate

of the jth marker in 3D space, cos

Ri=R(ai)

=

-sin

cYi cr,

i

, cos ff; 1 sin LY~

is the rotation matrix of the ith projection is the scale, Pi,

pi’ = [

Pi,

around the y-axis orthogonal

to the projection

(micrograph)

plane, si

1

are the original coordinates

of the jth marker on the ith projection,

and

1

dx, d

di= [

Y,

is the translation of the ith projection. Thus, given the marker coordinates pi the parameters ei, oi, si and di have to be found. Marker coordinates xj in 3D space are functions of marker coordinates in 2D and are not independent. In the algorithm proposed, they serve as a convenient model and as such are used in a second part of the method (Appendix A.21, which allows determination of the relative orientation of two single-tilt series. Eq. (1) can be solved only if the number of known parameters (2s~) is at least equal to the number of parameters to be determined (5~ - 2). Two of the orientation parameters are not independent since the rotation of the system around the y-axis is not fixed; thus one of the 0; angles can be arbitrarily set to zero. The same is true for the scale of the whole model, thus one of the scales si can be set to one. Thus 2s~ 2 5~ - 2, or s 2 3. In practice, there are errors in determination of markers’ locations on micrographs and, to minimize the overall error, the number of markers used should be larger. Overdetermination of the problem reduces the final error and allows rejection of incorrect markers. It is not required that all the markers should be present on all the micrographs. For the sake of simplicity we do not describe the necessary provisions for missing markers, but the modifications required are straightforward. The first step is the estimation of the locations of the markers’ coordinates in 3D space from their locations on the micrographs. We will show that the least-squares solution can be found using linear regression. Eq. (1) has the following expanded form: cos [

0

ei

0 1

-sin 0

ei

1[I [ Xj

;;=

IHH

sin (Y~ Plz i sin (Y~ cos (Yi ‘i + PY,

cos ffi

-

dx,

4,

(2)

P. Penczek et al. / Ultramicroscopy 60 (1995) 393-410

405

The right-hand side of Eq. (2) defines the positions of markers pi’ in 2D, taking into account orientation scale corrections. Eq. (2) can be written in explicit form as a system of equations:

and

xj cos Bi - z j sin Bi = p::,

(3)

# =pz, with i=

1, . . . . V; j=

1, . . . . S. The solution

of the second set of equations

is

(4) The first equation

has to be rewritten (assuming

+xi=

zj( _tane,)

that 8; f 7r/2 for i = 1, . . . , v):

p::’ cos

(5)

ei ’

which has the form of AjXi + Bj = yi’,

(6)

which shows that Eq. (5) forms s sets of linear equations, each consisting of v equations and each having 2 unknowns Aj ( = zj) and BJ ( = xj). Since Yj and Xi are known (they are functions of markers’ positions and tilt angles), the least-squares solution can be easily found for each set using linear regression. Together with the solution for JJ coordinates given by Eq. (4) they provide the coordinates of the markers in 3D space I’ describing a “3D marker model”. The availability of this “3D marker model” makes determination of corrections for all the parameters possible. Using the left side of Eq. (2) (or Eq. (3)) and the known tilt angles Bi we can calculate the positions of markers py’ in the projections. This is simply a projection of the marker in the direction of the given tilt image, as predicted by the current 3D marker model. These positions are used to find new orientation parameters (Y~,si and dj. The least-squares solutions for the orientation parameters are as follows: (1) Trunslufion ofprojection - determined by the minimum of the following functional: Ld(di)

Ip;:‘j-p;j-dil*,

= i j=

(7)

1

which has the solution

(8) (2) Scale - found using the so called symmetric error, which has an advantage that it does not depend on the choice of reference nor on the rotation angle between the two sets of points (for detailed discussion see Ref. [25]). The corresponding functional has the following form:

which can be written as L,( Si) = ;

$ I J=I

1Pq

+ si i j= I

Ip;‘1* - 2p;jTp;J.

(10)

406

P. Penczek et al. / Ulframicroscopy

Calculating

the derivative

60 (1995’) 393-410

over si and setting it to zero we obtain

which has a solution

1

[j=l

(12)

J.

j=l

(3) Rotation in plane - found by calculating

the minimum

of the functional

(13)

R( “j)p:jl*, (where R(oi) is rotation matrix) or, equivalently,

L,( cq) =

c

by finding

the maximum

of

p;‘R( q)p;j,

(14)

j=l

which can be expanded to L,( oi) = L

( p’:,p:: + p;ipz)

cos ai + i

j=l

Calculating

( pz/pz - piip:f) sin oi.

j=l

the derivative

over oi and setting it to zero we arrive at the following

solution:

Thus, we have shown that for parameters LY~,si and di we have a closed-form, least-squares solution. The only parameter for which such a solution cannot be found is the tilt angle 8;. The reason is that Eq. (5) is non-linear with respect to 8;. However, taking into account the fact that the tilt angles are known from the goniometer with a high accuracy, an approximate solution can be found. (4) Tilt angle - found using the following reasoning: projection of a point is achieved by the rotation of the system of coordinates with the third coordinate, p:, disregarded (see projection matrix f’;). Nevertheless, we can use the value of p: predicted by the 3D marker model (together with pz> to find the correction of the tilt angle. The rotation equation is as follows (we will omit the y-coordinate, since it does not play any role here):

(17) We know the coordinates of the marker in 3D (xj, z’> and its tilt angle 0;. Thus, using JZq. (17) we can calculate the predicted coordinates
( p;;)2+ ( p;;)’

=

( p:;)’+

( pqq2,

(18)

P. Penczek et al. / Ultramicroscopy

60 (1995) 393-410

407

and assuming that the point is on the same side of the plane as its corresponding point in the projection, we can obtain the sign of pi: and solve Eq. (18) in the form

(19) The correction for a current tilt angle fIi is found as an angle between pairs of points ( p::, p$ and ( p:(, pt,‘) using Eq. (16) with z-coordinates substituting x-coordinates. Using the solutions found the marker alignment algorithm is as follows: 1. using initial estimates for tilt angles and available marker coordinates, calculate corresponding 3D marker coordinates - 3D marker model (Eqs. (4) and (6)); 2. for an orientation parameter to be found, calculate predicted coordinates of markers in 2D from current 3D marker model and use them to find corrected value of orientation parameter (either of Eqs. (81, (12), (16) or procedure described in Section (4) >; 3. use the new orientation parameter to calculate the corrected 3D marker model; 4. repeat steps (ii> and (iii) for all the orientation parameters which have to be corrected; 5. repeat steps (ii>, (iii) and (iv> until the residual error of the 3D marker model becomes smaller than the preset value or cannot be reduced anymore. A.2. Orientation in 30 Based on the procedure described in Appendix A.l, two sets of marker coordinates in 3D space are known, one for each tilt series. Assuming that for each series the same markers were chosen, the orientation in terms of three Eulerian angles can be found. The coordinates of markers are x{ and xi, i = 1, . . . , s for the first and second set, respectively. The solution, optimal in a least-squares sense, is a rotation matrix R($,f?,q) (which is a function of three Eulerian angles) which maximizes the following expression: i j=

x:“R(‘//,8,p)x~.

(20)

1

Since Eq. (20) is highly non-linear and difficult to solve we will follow the solution suggested by Horn [25] instead. He used quatemion mathematics to find a closed-form least-squares solution (one that does not require iterations) of Eq. (20). A quatemion is a complex number with real part q0 and three imaginary parts q,, q, and q,: B=qO+iq,+jq,+~q,,

i*=

-1

j*=

ij=k, ji=

(21)

-k,

-1,

p=

-1;

jk= i,

ki=j;

kj= -i,

ik= -j.

(22)

Vectors in 3D space can be represented by purely imaginary quatemions, i.e., quatemions with q0 = 0. The product of two quatemions can be represented by matrix multiplication using the orthogonal 4 X 4 matrix: 40

-4,

-9)’

4,

40

-qz

qY

4z

40

-4y

4x

4,

-q*

qY -4, 40

H=Q2,

(23)

4.08

P. Penczek et a[./ Ultramicroscopy

60 (I9951 393-410

and 40

- 4-r

qv 4:

-qy

4:

-qv

-4,

40

4.x

4,

-4x

90

4,

-4:

K=U.i.

(24)

qo i

Rotation of a vector in 3D space can be represented

as a composite product of a purely imaginary

quatemion

K: K’ = $4 * = (Qi)G*

=d(a)

= (do)s.

(25)

Since 48 aTo=

0

0

0

2h,

(Y;+Y.:-Y:-q2

0

2(%&

+ 404,)

0

2(WX

- Yo4J

- 40%)

(d-d+Y;-YT) %?Zqy+404x)

0

GG% + 904,) 2(qd&- 404X) ’

(26)

(d-d-q,2+q;2)

and K is imaginary, therefore w’ has to be imaginary. In addition, since Q and -c) are orthonormal if C$is a unit quatemion (44 = 11, then lower-right-hand 3 X 3 submatrix of Eq. (26) must also be orthonormal and as such it is the rotation matrix. Thus, Eq. (25) defines a rotation. Using quatemions and the notation introduced we can write Eq. (20) in a different form:

c

(@$*).q,

(27)

j= I

where K{ and .?i are purely imaginary quatemions corresponding second set, respectively. Eq. (18) can be further modified

j=

1

to marker locations

in 3D in the first and

j= I

with the following

notation

introduced

0

_x{

-y;

4

0

-z:’

Y:

z:

-z; Y:

0

-Xi

4

-y:

x;

0

0

_xi

_y;

-z;

4

0

z:

Y:

-z:

0

x:

Y:

-Xi

0

z:

-Y:

(29)

P. Penczek et al. / Ultramicroscopy 60 (1995) 393-410

409

and (L

+ s,, + %> s,: - sz,

N=

I

SYZ- SZY (sYs,,-L)

s:.x- 5;

sx, + s,x

sx, - s,x

5, + sx:

L

- s,z

s,, - s,, s,* + sxz

SX” , + s,x

sy:+szy

(- sxx+ s,, - %) s,; + sz,

’

(-L-svy+szJ

where S,” = i j=

u;lJ:,

ll:,v:

=x;,y:,z:.

(30)

I

Thus, we showed that finding maximizes

the maximum

of Eq. (20) is equivalent

to finding

the unit quatemion

which

It can be shown (see Ref. [25]) that such a quatemion is the eigenvector corresponding to the most positive eigenvalue of the matrix N. The eigenvectors of a real symmetric matrix N can be found using standard numerical procedures [26]. The procedure of finding the rotation angles between two sets of markers x{ and xi, j = 1, . . . , s, in 3D space consists of the following steps: For each pair of coordinates all the possible cross-products are calculated and, using Eq. (301, the entries of the matrix N are calculated. The eigenvectors of matrix N arecalculated. From the eigenvector corresponding to the most positive eigenvalue, the 4 X 4 matrix is formed using Eq. (26), in which the lower-right-hand 3 X 3 submatrix is the required rotation matrix. Optionally, the corresponding three Eulerian angles can be retrieved from the entries of this submatrix using the procedure described in Ref. [27]. The rotation found in step (3) is applied to the second set of marker coordinates and the residual error between each pair of markers from the respective sets is calculated to estimate the quality of fit.

References ill CA. Mannella, M. Marko, P.A. Penczek, D. Barnard and J. Frank, Microsc. Res. Tech. 27 (1994) 278. l21 B. McEwen, J.T. Arena, J. Frank and CL. Rieder, J. Cell Biol. 120 (1993) 301. 131 R.A. Horowitz, D.A. Agard, J.W. Sedat and CL. Woodcock, J. Cell Biol. 125 (1994) 1. [41 A.L. Olins, D.E. Olins, H.E. Levy, M.B. Shah and D.P. Bazett-Jones, Chromosoma 102 (1993) 137. [51 A.S. Belmont, Y. Zhai and A. Thilenius, J. Cell Biol. 123 (1993) 1671. l61 K. Schmeckel, J. Wahraman, U. Skoglund and B. Daneholt, Chromosoma 102 (1993) 669. 171 M.S. Ladinsky, J.R. Kremer, P.S. Furcinitti, J.R. McIntosh and K.E. Howell, J. Cell Biol. 127 (1994) 29. ls1 M. Marko, P. Penczek, R. Feuerbach and J. Frank, in: Proc. 13th Int. Congr. on Electron Microscopy, Paris, 1994, Vol. 1 (Les Editions de Physique, Les Ulis, 1994) p. 505. [91 M. Radetmacher, J. Electron Microsc. Tech. 9 (1988) 359. [lOI W. Hoppe- and R. Hegerl, in: Computer Processing of Electron Microscope Images, Ed. P.W. Hawkes (Springer, Berlin, 1980). Imaging with the Transmission Microscope, Ed. J. Frank (Plenum, New [ill J.M. Carazo, in: Electron Tomography: Three-Dimensional York, 1992). [I21 D.F. Parsons, M. Marko and M.V. King, J. Microscopy 118 (1980) 127. 24 (1988) 7. [I31 P.K. Luther, MC. Lawrence and R.A. Crowther, Ultramicroscopy [141 L.A. Hindus, Advanced Imaging 7 (1992) 52. in: Electron Tomography: Three-Dimensional Imaging with the Transmission Microscope, Ed. J. Frank (Plenum, [151 M. Radennacher, New York, 1992).

410 1161

P. Penczek et al./

Ultramicroscopy

60 (1995)

393-410

P.A. Penczek, M. Radetmacher and J. Frank, Ultramicroscopy 40 (1992) 33. [17] W.O. Saxton, W. Baumeister and M. Hahn, Ultramicroscopy 13 (1984) 57. [18] M.A. Levy, A.L. Olins and D.E. Olins, J. Microscopy 165 (1992) 325. [ 191 K. Dierksen, D. Typke, R. Hegerl and W. Baumeister, Ultramicroscopy 49 (1993) 109. [20] K. Dierksen, D. Typke, R. Hegerl and J. Walz, Biophys. J. 68 (1995) 1416. [21] J. Frank, B. Shimkin and H. Dowse, Ultramicroscopy 6 (1981) 343. [22] MC. Lawrence, Proc. Electron Microsc. See. Afr. 13 (1983) 19. [23] M.C. Lawrence, in: Electron Tomography: Three-Dimensional Imaging with the Transmission Microscope, Ed. J. Frank (Plenum, New York, 1992). [24] J. Berriman, R.K. Bryan, R. Freeman and K.R. Leonard, Ultramicroscopy 13 (1984) 351. [25] B.K.P. Horn, J. Opt. Sot. Am. A 4 (1987) 629. [26] W.H. Press, S.A. Teukolsky, W.T. Vetterling and B.P. Flannery, Numerical Recipes (Cambridge University Press, Cambridge, 1992). [27] P.A. Penczek, R.A. Grassucci and J. Frank, Ultramicroscopy 53 (1994) 251.