key: cord-0935806-3u81ntle authors: Acar, Erman; Peltonen, Sari; Ruotsalainen, Ulla title: Adaptive multiresolution method for MAP reconstruction in electron tomography date: 2016-08-06 journal: Ultramicroscopy DOI: 10.1016/j.ultramic.2016.08.002 sha: 2e2b5f57b8dad06f44a672757c71fdc37802420f doc_id: 935806 cord_uid: 3u81ntle 3D image reconstruction with electron tomography holds problems due to the severely limited range of projection angles and low signal to noise ratio of the acquired projection images. The maximum a posteriori (MAP) reconstruction methods have been successful in compensating for the missing information and suppressing noise with their intrinsic regularization techniques. There are two major problems in MAP reconstruction methods: (1) selection of the regularization parameter that controls the balance between the data fidelity and the prior information, and (2) long computation time. One aim of this study is to provide an adaptive solution to the regularization parameter selection problem without having additional knowledge about the imaging environment and the sample. The other aim is to realize the reconstruction using sequences of resolution levels to shorten the computation time. The reconstructions were analyzed in terms of accuracy and computational efficiency using a simulated biological phantom and publically available experimental datasets of electron tomography. The numerical and visual evaluations of the experiments show that the adaptive multiresolution method can provide more accurate results than the weighted back projection (WBP), simultaneous iterative reconstruction technique (SIRT), and sequential MAP expectation maximization (sMAPEM) method. The method is superior to sMAPEM also in terms of computation time and usability since it can reconstruct 3D images significantly faster without requiring any parameter to be set by the user. Electron tomography (ET) is a widely used powerful technique in the field of biology revealing the interior structure of the biological samples in 3D at nanometer scale. This is achieved by collecting 2D projections of the specimen from different viewing angles. In its common application, the sample is rotated around a single axis with fixed intervals while its projections are acquired by transmission electron microscope. Then, these projections are used to reconstruct the biological sample in 3D. The accuracy and the resolution of the reconstruction are limited by two major factors: (1) due to the mechanical and physical limitations, the projections typically collected only from 7 60°to 70°tilt angle range are used in the reconstruction (missing wedge problem). (2) The signal to noise ratio (SNR) of the projection data is low since the total amount of applied electrons is limited in order to prevent radiation damage to the biological sample imaged. The maximum a posteriori (MAP) reconstruction methods are well suited for reconstruction from such incomplete and noisy projection data. They are able to use a priori information to compensate for the imperfection and incompleteness of the measurements. Any structural or statistical information about the sample or the imaging conditions can be utilized to improve the reconstruction process. These methods try to satisfy the fidelity of the reconstruction to the measurements and utilize the prior information at the same time. The balance between the data fidelity and the prior information is controlled by the regularization parameter. The value of this parameter highly affects the visual quality of the reconstructed image. Its optimum value depends on several factors, such as the noise level, sample variation, resolution, and validity of the prior information. The computation of the regularization parameter in MAP estimation is usually based on the SNR analysis of the data, structural information about the objects reconstructed, or predefined set of parameters aiming to provide a solution as general as possible. One of these methods, the L-curve method, analyzes the images reconstructed for different values of regularization parameter [1] . For each reconstruction, the relationship between the squared norm of the solution and the squared norm of the residual error is analyzed to find the optimum regularization parameter. A similar approach, S-curve method, was implemented for wavelet-based reconstruction [2] and total variation regularization [3] making the assumption that the reconstructed data is sparse in wavelet domain and gradient domain, respectively. Both L-curve and S-curve methods require multiple reconstructions of the data with different regularization parameters, which heavily increases the computational cost of the overall image reconstruction. In another method, Wen and Chan used discrepancy principle for total variation regularization in the image restoration context [4] . Their method requires prior information about the SNR of the data, which is not always straightforward to obtain in tomographic imaging. In this study, one of our aims was to develop an adaptive regularization method, which does not require additional analysis, structural a priori information about the data nor the imaging system, or the user manipulation according to the application purpose. The other aim of this paper was to reduce the computation time, which is realized by employing multiresolution reconstruction in a sequential scheme. Multiresolution reconstruction grids [5] and multiresolution detector space [6] was used previously to improve the convergence rate of maximum likelihood expectation maximization (ML-EM) method in positron emission tomography (PET). It was also utilized in scanning transmission electron microscopy (STEM) [7] and discrete tomography that reconstructs images composed of small number of intensity levels representing only a few different materials [8] . In this paper, we use the multiresolution method in electron microscopy for the reconstruction of biological specimens with maximum a posteriori probability expectation maximization (MAPEM) method. Multiresolution approach introduced in MAPEM reconstruction enables the weight of the regularization to vary from strong to weak throughout the reconstruction. One reason is that the image size is increased during the reconstruction while the regularization filter size is fixed. The other reason is the number of binned projection data and so its averaging effect decreases throughout different resolution stages. By decreasing noise with these two mechanisms, the reconstruction converges smoothly to a coarse estimate of the reconstructed image at the first resolution stages. Then the number of binned pixels is decreased step by step to improve the resolution of the estimates. Since each step is initialized with a spatially coarse estimate of the image, they give better results in fewer iterations compared to initialization with an image of random or constant value. The computations take shorter time also because the calculations are faster in the smaller reconstruction grids used in the lower resolution stages. We developed the adaptive regularization method first as a parameter-free reconstruction process and then introduced the multiresolution approach to obtain the final adaptive multiresolution reconstruction method. An accurate model for the formation of TEM images would be based on quantum and relativistic aspects of the illumination of the specimen, electron-specimen interaction, microscopy optics, and the detection of the electrons. However, this would be a computationally unfeasible model. Therefore, we apply some simplifications and statistical modeling to obtain a feasible image formation model. A transmission electron microscope can be operated in different imaging modes depending on the illumination of the specimen, the contrast mechanism based on the interaction of the electrons with the specimen, and the detector type. In the dark field imaging mode, the detection of scattered electrons plays the major role in image formation, whereas it is the detection of transmitted electrons in bright-field imaging mode. We consider the bright field mode, used conventionally for ET in life sciences, in this study. In this imaging mode, a parallel beam of electrons illuminates the specimen, and the image is formed in the back focal plane of the objective lens by the electrons that pass through the aperture. For the image contrast mechanism, we only consider the amplitude contrast generated due to the electrons that are blocked from reaching the detector and we ignore the phase shift created by the electron-specimen interaction. This model is sufficient for medium-resolution contrast (beyond 2-3 nm) imaging however phase contrast should also be considered at higher resolutions where the features to be detected are smaller than the coherence length of the electron [9] . We also assume that successive electrons do not interact with each other so that they can be treated independently. This assumption holds for TEM imaging of biological specimens since the specimen is thinner than the mean distance between two successive electrons [10] . Under these assumptions the formation of the image can be expressed as: is the data acquired for the tilt angle θ, C is the data that would be acquired if there was no specimen, PSF represents the point spread function of the imaging system, * is the convolution operation, and ( ) Here, v b is proportional to the probability that an incident electron traveling along θ is scattered inside the bth voxel of the discretized specimen volume. The contribution of the bth voxel to the dth projection data, p d , is represented by a db . C and PSF in (1) can be determined by separate calibration experiments and they can be excluded from the equation by subtraction and deconvolution. After this basic preprocessing of the raw data, the problem reduces to estimate the parameters of distribution describing v b 's in (2) . In parallel beam projection geometry, the reconstruction space can be divided into slices orthogonal to the tilt axis and the slices can be reconstructed independently. Then the slices can be concatenated to reconstruct the final 3D volume. Therefore, the method will be described here for single slice of the volume to be reconstructed. The image formation geometry is presented in Fig. 1 . The parameter estimation can be realized for each slice independently by MAPEM method which maximizes a posteriori probability distribution using the Baye's rule as: ) P Data Img is the likelihood of the measured data to be observed given the image to be reconstructed. ( ) P Img and ( ) P Data are the a priori probability distributions. ( ) P Data can be assumed to be constant and excluded from the expression. In ET, very low dose of electrons is applied to a very thin sample in order to prevent radiation damage to the biological sample. Considering the characteristics of the signal and the noise, Poisson distribution can be assumed to represent the measurements. This assumption can be extended to model the number of electrons blocked by each pixel of the image. The statistical chain from generation of the electrons, via the scattering events inside the specimen, to the measurement detector, is assumed to be a cascaded Poisson process. This assumption is useful to satisfy the positivity constraint throughout the iterations. It was also utilized in the development of the adaptive regularization method described below (in Section 2.2). Using Poisson distribution and (2), the likelihood of observing the projection data, p, conditioned on the mean values λ ( | ( ) P Data Img ) can be expressed as: where λ b is the parameter of the Poisson distribution describing v b . By taking the logarithm and derivative of this expression with respect to λ b and equating it to zero, the maximum likelihood can be obtained as: Since v a b db is unavailable, it will be replaced by its conditional expected value at each iteration of the expectation maximization algorithm. The expected value can be calculated as: where ′ b is the pixel index of the image to be reconstructed. Using the expected value of v a b db in (4), the maximum likelihood estimate of λ can be expressed iteratively as follows. In addition to the likelihood term, C b L k , MAPEM reconstruction involves a regularization term, C b P k , that represents ( ) P Img . The overall update equation of the MAPEM is For ( ) P Img , local similarity in the small neighborhood of each pixel of the reconstructed image is a common choice. With Median Root Prior (MRP) [11] , the regularization term is expressed as: where () med operator represents the median filtering and β is the regularization parameter which controls the strength of this penalization filter. Sequential maximum a posteriori expectation maximization (sMAPEM), was introduced to ET recently to compensate for the missing wedge effects [12, 13] . The method assumes Poisson distribution to model the image to be reconstructed and median filter (median root prior, [11] ) to regularize the iterations. Instead of using one constant regularization parameter value, the method used a sequence of these in the decreasing order. At each regularization stage, predefined number of iterations was performed with constant regularization parameter value. The initial stages with high values of regularization parameters provided a rough but robust estimate of the image for the following stages by highly suppressing the noise. The regularization parameter value was decreased to enhance the resolution and contrast of the final reconstructed image. This sequential method was consistently superior to the conventional reconstruction methods, weighted back-projection (WBP) and a simultaneous iterative reconstruction technique (SIRT), in terms of accuracy and noise suppression due to its capability of filling the missing wedge [12] . However, sMA-PEM has some parameters, which affect the final resolution, the final contrast, and the visual impression of the reconstruction result. The number of regularization stages, the number of iterations 0º Fig. 1 . The projection geometry. The volume is divided into x-z slices which are reconstructed independently. v b is proportional to the probability that an incident electron is scattered inside the bth pixel of the x-z slice image. It is modeled by a Poisson distribution with mean λ b . p d is the dth projection data. per each stage, and the regularization parameters need to be set properly by the user according to the data and the application. Moreover, the method requires a long computation time. MAPEM methods maximize the likelihood of observing the measured projections subject to the prior knowledge. During this maximization, noise in the projection data makes the pixels with high variance of intensity values vary around their expected values more than the ones with low variance. It would be wise to adjust the regularization parameter pixel-wise considering the amount of local variance. This way, the prior information can be used more strongly for the pixels with high variance of intensity values than the ones with low variance of intensity values. The intensity values of the pixels at each iteration correspond to the mean values of the Poisson distributions. For Poisson distribution, the variance of the distribution is equal to its mean value. Therefore, the intensity values at each iteration, λ b k , can be used as a measure of variance and so for the regularization strength. For the positivity constraint of the pixels, the regularization parameter, β, needs to be in the range [0,1]. Therefore, the intensity values were normalized considering this range at each iteration and used as the pixel-wise regularization parameter in this study. Hence, the adaptive regularization term can be expressed as: The likelihood term, C b L k , and the overall update equation are the same as in (6) and (7), respectively. Note that, we assign a different regularization parameter to each pixel at each iteration instead of using a constant regularization parameter for the whole image. The strength of the penalization depends on the regularization parameter value and the filter size relative to the reconstructed image size. sMAPEM method used a sequence of regularization parameter values in the decreasing order to weaken the penalization filter sequentially [12] . In this way, noise was highly suppressed in the reconstruction in the first regularization stages to have coarse but robust estimate of the structure. In the following stages, the regularization parameter value was decreased to enhance the resolution and contrast of the estimate. In addition to changing the regularization parameter value adaptively to control the noise penalization, adaptive multiresolution MAPEM (amMA-PEM) reconstruction method changes the reconstruction image size while keeping the filter size constant. In our experiments, we fixed the window size of the filter to its smallest value (3  3) in order to minimize the loss of resolution due to filtering (The reader may refer to [11] for the effect of different filter sizes on MRP regularization). The same size of filter affects the small size image more strongly than the large one. With amMAPEM method, the strength of the penalization filter is weakened by increasing the image size stage by stage sequentially. The advantage of this approach is its computational efficiency. The flow diagram of the amMAPEM method is given in Fig. 2 for single x-z slice of the reconstructed volume. The method is initialized with a uniform image of size 4  4 so that 3  3 median filter can fit in properly. At each resolution stage, adaptive regularization is applied using the binned projection data and rescaled reconstruction grid. The result of each regularization stage is scaled up by a factor of 2 using bilinear interpolation. The scaled image initializes the following regularization stage. The binning factor of the projection data is also decreased throughout the stages. With this method, the averaging effect of binning operation decreases the noise contamination together with the median filter employed inside the iterations of adaptive regularization. Decreasing the number of binned pixels step by step improves the resolution of the reconstruction. The final image is obtained at the last stage using the result of previous stages and the original projection data P. The method is implemented in MATLAB (MathWorks Inc., MA, USA) and the code is available together with the data related to this study at http://www.cs.tut.fi/sgn/m2obsi/m2obsiWWW/de mos/amMAPEM/amMAPEM.html. The iterations for each resolution stage of amMAPEM were stopped when the difference between two consecutive images throughout the iterations was insignificant. The difference was measured by normalized mean squared error (NMSE) and the level of significance was determined to be − 10 7 experimentally. The NMSE values are calculated as: A 256  256  64 nm 3 simulated numerical phantom was generated in MATLAB for the evaluation of the methods. It is composed of about 25 spherical objects located additively in the volume. 25 hollow spheres were also added to the phantom to simulate vesicles. The intensity levels and diameters of the objects were determined randomly following a uniform distribution in the range [0 1] and [19, 42] nm, respectively. The projection images were calculated with 1°angular step in the range 760°. The images were generated at a pixel size of 1 nm  1 nm. Two types of noise were added to the projection data (1) Poisson noise to simulate the randomness in the counting process of transmitted electrons and (2) Gaussian noise to simulate the variations at the detector level. Three noise levels were simulated to have SNRs (ratio of the signal and the noise variances) of 50, 10, and 1. The zero tilt images of the dataset are presented in Fig. 3 . The experimental dataset consists of a freeze-substituted Ver-oE6 cultured cell infected with the SARS-CoV (The Cell Centered Database, project P2005, microscopy product 6021, [14, 15] ). It includes 131 TEM images acquired with a 1°angular step at the range [ À 65°, þ65°] using single tilt projection geometry. The images were collected at 80 kV with a Philips CM-10 transmission electron microscope (Philips, The Netherlands). Spherical gold particles with diameter 10 nm were applied to the surface of the specimen to serve as fiducial markers to align the tilt series. The original images with size 2048  2048 were cropped to 1024  1024 pixels for computational simplicity. The pixel size was 1.2 nm. The zero tilt image of the dataset is presented in Fig. 3d . The experimental dataset is comprised of a single axis tilt series of chlamydia trachomatis type III secretion systems in contact with a HeLa cell [16, 17] . It contains 36 cryoEM images collected in 3°i ncrements at the range [ À 60°, þ45°]. The dataset also includes spherical gold particles with diameter 10 nm. The alignment and CTF correction of the projection images were realized in IMOD [18] . The raw projection image size was 3708  3838 pixels. 2  2 binning and cropping were applied to the images after the alignment and CTF correction process. The resulting image size was 1369  1448 and the pixel size was 1.08 nm. The zero tilt image of the dataset is presented in Fig. 3e . In order to have a general assessment of the overall image quality, the mean squared error (MSE) between the reconstruction result and the ground truth was calculated as: N is the number of voxels, I gt is the ground truth and I is the reconstructed volume. Fourier shell correlation (FSC) is the normalized cross-correlation coefficient between two 3D volumes calculated over shells in Fourier space [19] . It can be calculated using the ground truth and the reconstructed volume as: where F gt and F are the Fourier transforms of the ground truth and the reconstructed volume, respectively. r i is the ith voxel at radius r. FSC expresses the correlation between two volumes with respect to spatial frequency. The cutoff frequency providing correlation higher than a certain threshold value is used as a measure of resolution in the field of electron microscopy. In order to evaluate the contribution of adaptive multiresolution method to the reconstruction process, the evaluation results were compared with sMAPEM, and the conventional methods, WBP and SIRT [20] . sMAPEM reconstructions were performed using the recommended parameters. In WBP reconstructions, a ramp filter and a Hamming filter with cutoff frequency 0.5 cycles/pixel were used. 50 iterations were used for the SIRT reconstructions. Positivity constraint was applied for WBP and SIRT reconstructions. All the methods were implemented in MATLAB. The MSE values calculated between the reconstruction results and the ground truth are presented in Table 1 for different noise levels. The MSE value of amMAPEM was lower than the WBP, SIRT and sMAPEM values on average. The values for different noise levels show that amMAPEM is more robust against noise than sMAPEM. The increasing noise results in bigger increase in MSE for sMAPEM than for amMAPEM. We also calculated Fourier shell correlation (FSC) between the reconstructed images and the ground truth to evaluate the images for their capability of resolving the details. The FSC curves are shown in Fig. 4 . The correlation values get close to 1 as the reconstruction result gets close to the ground truth. It is seen in the figure that the correlation values with respect to spatial frequencies are clearly higher for sMAPEM and amMAPEM than for WBP, SIRT. The correlation with the ground truth decreases with the increasing noise. The decrease is bigger for sMAPEM compared to amMAPEM. This shows amMAPEM can handle noise variations better than sMAPEM. The visual assessment of the reconstruction results is in line with the numerical evaluations. The resulting images of the reconstructions for the highest (SNR 1) and lowest (SNR 50) noise levels are given in Fig. 5 . The missing wedge artifacts are clearly visible in WBP and SIRT reconstruction results. The spherical objects are elongated in z direction and the elongated objects create artificial objects in the x-y slices. It is also seen that hollow spheres are deformed and look open in the in the WBP and SIRT reconstruction results because of the missing wedge. However, these missing wedge effects are compensated better in sMAPEM and amMAPEM reconstructions. Fig. 5 also shows that the noise suppression capability of the amMAPEM method with the adaptively selected parameters is better than that of the sMAPEM method. The effect of noise is visible inside the smooth regions of the spherical objects. The results also show that the amMAPEM method can preserve the edges while suppressing the noise. It is clear in the images that the smooth regions are smoother with amMAPEM than with sMAPEM while the boundaries are still sharp. Fig. 6 shows the line profiles of a circular object in z direction. WBP and SIRT profiles are wider than the ground truth profile while the sMAPEM and amMAPEM profiles show edges as sharp as the ground truth. Both sMAPEM and amMAPEM can keep the object boundaries quite close to the ground truth. The major difference between the profiles of these two methods is seen in the plateau with the high intensity values. amMAPEM profile is flatter than sMAPEM profile in this region for all noise levels due to the adaptive regularization. The slice images taken from the reconstruction results for different methods are analyzed in Fig. 7 . The zoomed-up views of the gold particle show the elongation effect in z direction due to the missing wedge. The spherical shape of the gold particle is highly distorted in WBP and SIRT reconstruction results. However, it is still close to its original spherical shape in the sMAPEM and am-MAPEM reconstructions. In Fig. 8 , the gold particle images are analyzed in more detail using line profiles taken through the center of the particle shown in Fig. 7 . The profiles taken in the x and z directions show that the edges of the sMAPEM and am-MAPEM profiles are sharper than the WBP and SIRT profiles. The plateau with the high intensity values is supposed to be flat because of the homogeneity of the gold particle. This region is flatter for amMAPEM compared with sMAPEM. The gold particles in the reconstructed volumes were also analyzed quantitatively. For that purpose, randomly selected 5 gold particles were isolated and compared in terms of FSC with the ground truth generated using the known size of the gold particles. The x-z and x-y projection images for three of these gold particles are presented in Fig. 9 . The average of 5 FSC curves is presented in the figure. The images and the FSC curve show that sMAPEM and amMAPEM methods provide more accurate results with less elongation artifacts compared to WBP and SIRT. The slice images taken from the reconstruction results are presented in Fig. 10 . The WBP images are noisier than the other images. The boundary of the circular object is preserved better in amMAPEM image compared to the other results. The quantitative analysis of the reconstructed volumes was performed using the gold particles as it was done for the CCDB-P2005 dataset. The x-z and x-y projection images for three of the isolated gold particles are presented in Fig. 11 . The average of FSC curves calculated using Ground Truth xz slices xy slices xz xy The datasets were divided into x-z slices and the slices were reconstructed independently. WBP and SIRT slices were reconstructed one by one using a desktop PC. For amMAPEM and sMAPEM, a grid of computers with various configurations was utilized to reconstruct the slices in parallel. In order to compare the computational efficiencies of the reconstruction methods, we measured single slice computation times for each method. The computations were performed using MATLAB on an Intel(R) Core (TM) i5@3.20 GHz desktop computer with Windows 64 bit operating system. The computation times are presented in Table 2 . The computation time for amMAPEM increases with the increasing noise. The computation times of the other methods are same for all noise levels since they use fixed parameter settings. amMAPEM was significantly faster than the sMAPEM method. Although the WBP and SIRT computation times were shorter than those of sMAPEM and amMAPEM, the visual quality and accuracy of these reconstructions were much worse than sMAPEM and amMAPEM. The present research demonstrates that the proposed amMA-PEM method can reconstruct 3D ET images more accurately than the conventional reconstruction methods, WBP and SIRT, significantly suppressing the missing wedge artifacts. It also provides images with better quality than with sMAPEM in shorter time without requiring special effort for parameter setting. The proposed adaptive regularization method extracts all necessary information simultaneously from the reconstructed image throughout the iterations. Instead of using a single regularization parameter for the whole image, it updates the regularization parameter pixel-wise according to the local noise contamination. This way, we calculate the regularization strength adaptively. The major advantage of sMAP-EM and amMAPEM methods is their gap filling ability. The methods provide a compensation for the missing angular gap in the frequency domain by using the prior information of local smoothness. The gap filling ability arises from the MRP regularization. The median filter used throughout the iterations acts as an interpolator in the frequency domain. The information source for the interpolation is the relation of the pixel values in the local neighborhood. The filter preserves the edges in the spatial domain while using this information. The reader may refer to the previous studies for the detailed analysis of gap filling capability of sMAPEM in comparison with the conventional reconstruction methods [12, 13] . One of the major drawbacks of sMAPEM is that it requires a large number of iterations for gap filling. We overcome this computation time problem with am-MAPEM method by employing multiresolution reconstruction. The amMAPEM method uses the intensity values, which represent the mean of the Poisson model used for each pixel, as a measure of noise contamination since the variance of a Poisson distribution equals its mean. With Poisson model, it was possible to characterize the distribution by a single positive parameter. This provided positivity constraint to be satisfied implicitly for the reconstruction in addition to simplicity of the solution. The usability of the iterative reconstruction methods in electron tomography highly depend on their computation time. The opportunities provided by the modern computing technology (e.g. GPU [21] , SIMD extensions of the modern processors [22] ) can be utilized to improve the performance of the methods. We used the computer grid Techila (Techila Technologies Ltd., Tampere, Finland) for the reconstruction of the experimental data in this study. However, methodological improvements still play an important role in accelerating the image reconstruction besides these improvements in the software implementation level. Ordered subsets (OS) is one method to decrease the computation time in tomographic image reconstruction [23] . The method divides the projection data into predefined number of subsets and updates the reconstructed image for each subset during the iterations. Since there is a large block of information gap and SNR is low in ET, updating the image using even sparser subsets of projection data can exaggerate noise throughout the iterations. Therefore, we utilized the multiresolution approach instead of OS to reduce the computation time. It was shown with the phantom and experimental datasets that the computation time is much shorter for amMAPEM than for sMAPEM. The multiresolution reconstruction also contributed to the visual quality of the reconstruction results in addition to computation time acceleration. Its averaging effect and usage together with the median filter suppressed noise throughout the iterations. Each regularization 11 . The gold particle based analysis of the experimental data EMPIAR-10048 reconstruction results. The x-z and x-y projection images of three isolated gold particles are presented on the left. The average FSC curves calculated between the reconstruction results and the gold particle ground truth is presented on the right. sMAPEM and amMAPEM methods provide more accurate reconstruction results with less elongation artifacts compared to WBP and SIRT. sMAPEM images look noisier than amMAPEM images. stage (except the first one) is initialized with a coarser estimate of the image. Starting with this coarser estimate of the image to be reconstructed instead of a random or uniform one improved the reconstruction result. The MAP reconstruction method presented in this paper provides a successful solution to compensate for the missing wedge artifacts and noise in electron tomography. Since the method does not require any parameter setting and it has improved computational efficiency, it can be widely used in the field of electron microscopy. Rank-deficient and discrete ill-posed problems: numerical aspects of linear inversion Sparse tomography Multi-resolution parameter choice method for total variation regularized tomography Parameter selection for total-variation-based image restoration using discrepancy principle A multigrid expectation maximization reconstruction algorithm for positron emission tomography Multiresolution expectation maximization reconstruction algorithm for positron emission tomography using wavelet processing A model based iterative reconstruction algorithm for high angle annular dark field-scanning transmission electron microscope (HAADF-STEM) tomography A multiresolution approach to discrete tomography using DART Handbook of Mathematical Methods in Imaging Transmission Electron Microscopy and Diffractometry of Materials Generalization of median root prior reconstruction Compensation of missing wedge effects with sequential statistical reconstruction in electron tomography A Bayesian approach for suppression of limited angular sampling artifacts in single particle 3D reconstruction A cell-centered database for electron tomographic data SARS-coronavirus replication is supported by a reticulovesicular network of modified endoplasmic reticulum EMPIAR: a public archive for raw electron microscopy image data Structure of a bacterial type III secretion system in contact with a host membrane in situ Computer visualization of threedimensional image data using IMOD Fourier shell correlation threshold criteria Fundamentals of three-dimensional reconstruction from projections High-performance iterative electron tomography reconstruction with long-object compensation using graphics processing units (GPUs) Vectorization with SIMD extensions speeds up reconstruction in electron tomography Accelerated image reconstruction using ordered subsets of projection data This work was supported by Tekes FiDiPro (1913/31/2012). The experimental dataset was provided by The Cell Centered Database which is supported by NIH Grants from NCRR RR04050, RR RR08605 and the Human Brain Project DA016602 from the National Institute on Drug Abuse, the National Institute of Biomedical Imaging and Bioengineering and the National Institute of Mental Health.