Abstract
Electron ptychography has seen a recent surge of interest for phase sensitive imaging at atomic or nearatomic resolution. However, applications are so far mainly limited to radiationhard samples, because the required doses are too high for imaging biological samples at high resolution. We propose the use of nonconvex Bayesian optimization to overcome this problem, and show via numerical simulations that the dose required for successful reconstruction can be reduced by two orders of magnitude compared to previous experiments. As an important application we suggest to use this method for imaging single biological macromolecules at cryogenic temperatures and demonstrate 2D singleparticle reconstructions from simulated data with a resolution up to 5.4 Å at a dose of 20e ^{−}/Å^{2}. When averaging over only 30 lowdose datasets, a 2D resolution around 3.5 Å is possible for macromolecular complexes even below 100 kDa. With its independence from the microscope transfer function, direct recovery of phase contrast, and better scaling of signaltonoise ratio, lowdose cryo electron ptychography may become a promising alternative to Zernike phasecontrast microscopy.
Introduction
The advent of direct electron detectors has led to a resolution revolution in the field of cryo electron microscopy in the last few years. The technique is now producing threedimensional atomic potential maps of biological macromolecules of a few 100 kDa or lower with a resolution better than 3.5 Å^{1,2,3}, such that individual amino acid sidechains can be resolved. An important role in this revolution play new image processing algorithms based on a Bayesian approach, which infer important parameters without user intervention^{4}. Also the correction of beaminduced motion has become possible mostly due to the new generation of detectors^{5, 6}. However, several challenges remain to be overcome in order to routinely reach 3 Å resolution also for small complexes^{7, 8}: Firstly, beaminduced specimen charging and subsequent motion currently still render the high resolution information of the first few frames of a high repetition rate movie recorded with a direct electron detector unusable^{9}, because the motion is too fast to efficiently correct for it. Secondly, the Detective Quantum Efficiency (DQE) of detectors is still imperfect at high spatial frequencies^{8, 10}. Thirdly, the contrast of single images can still be improved to enable reconstructions with fewer particles and increase the throughput^{10}.
The last of these challenges has recently been addressed with a new phase plate model^{11, 12}, which is comparatively simple to use and provides excellent contrast at low spatial frequencies. In addition to this hardwarebased approach to achieve linear phase contrast in the measured amplitudes, discovered by Zernike in the 1930s^{13}, it is also possible to algorithmically retrieve the phase information from a set of coherent diffraction measurements. One such technique, commonly known as ptychography or scanning coherent diffractive microscopy ^{14}, is becoming increasingly popular in the field of materials science due to experimental robustness and the possibility to obtain quantitative phase contrast over an essentially unlimited field of view^{15, 16}. The use of ptychography for imaging radiation sensitive samples with electrons at high resolution is however precluded so far by its high dose requirements.
Here, we show how the use of nonconvex Bayesian optimization to solve the ptychographic phase retrieval problem fulfills the dose requirements for imaging biological macromolecules and makes it possible to obtain 2D images from single particles with subnanometer resolution. After a short introduction into the technique, we will also mention how ptychography offers improvements for the other two challenges discussed above.
Despite being initially proposed as a solution to the phase problem for electrons^{17, 18}, ptychography has seen its biggest success in Xray imaging, due to the less stringent sample requirements and the experimental need for lensless imaging techniques. Recent developments include the introduction of iterative algorithms to enable the reconstruction of datasets collected with an outof focus probe^{19, 20}, which decreases the memory requirements of the method dramatically. The algorithms also have the capability for the correction of experimental difficulties such as unknown scan positions^{21,22,23}, partial coherence^{24}, probe movement during exposure^{25, 26}, intensity fluctuations during the scan^{24, 27} and reconstruction of background noise^{15, 27}.
In recent years, some of these advances have been applied in the context of electron microscopy and yielded atomic resolution reconstructions of lowatomic number materials^{28,29,30} and quantitative phase information^{15}.
Figure 1 shows the experimental setup for an outof focus ptychography experiment. A ptychographic dataset is collected by scanning a spatially confined, coherent beam, subsequently called ‘probe’, over the specimen and recording farfield diffraction patterns at a series of positions such that the illuminated regions of neighboring positions overlap. The diffractionlimited resolution r _{d} (halfperiod) of the final image is given by \({r}_{{\rm{d}}}=\frac{\lambda \cdot {\rm{\Delta }}z}{{N}_{{\rm{pix}}}\cdot {d}_{{\rm{pix}}}}\), where N _{pix} is the number of detector pixels, d _{pix} is the detector pixel size and λ is the deBrogliewavelength of the electrons. Given the set of positions \({\vec{r}}_{{\rm{1}}\ldots {\rm{K}}}\) and a realistic forward model for the formation of the corresponding diffraction patterns \({I}_{(1\ldots K)}\), the complexvalued transmission function \(T(\vec{r})\), which describes the atomic properties of the specimen^{31}, can then be retrieved by solving a nonconvex inverse problem. The electron dose used for successful reconstruction has exceeded 1 × 10^{3} e ^{−}/Å^{2} so far, limiting the usability of ptychography to radiationhard specimens. Table 1 lists recently published electron ptychography experiments and the used average electron doses. The lowest dose was reported in ref. 32, which used an estimated 3.33 × 10^{3} e ^{−}/Å^{2} at a resolution of d_{r} = 58.4 pm, resulting in a dose of 1.1 × 10^{3} e ^{−}/pixel and achieving a line resolution of 2.3 Å, demonstrated by resolving the lattice spacing of gold nanoparticles. We will demonstrate via simulations that it is possible to reduce this dose by a factor of 100, thus reaching the dose range allowed for imaging biological macromolecules.
The problem of beaminduced sample movement has already been addressed before the development of fast direct detectors. Scanning with small spots of several 10 nm in size over a vitrified sample has shown to reduce beam induced specimen movement^{33,34,35} in realspace imaging and it has been noted that the remaining movement may be due to radiation damage, not sample charging^{33}. Ptychography naturally operates with a confined beam, thus minimizing the area where charge can build up, such that the movement should be reduced compared to the illumination of large areas in cryoEM. The sampling requirements given by the experimental setup allow to operate the detector with large effective pixel sizes, such that the DQE and MTF are at nearperfect values and can be neglected in the reconstruction. The overhead resulting from the need to take multiple exposures can therefore be greatly reduced by using fast detectors with few pixels.
Results
Image formation of cryoTEM and ptychography
We perform multislice simulations of three different biological macromolecules with molecular weights ranging from 64 kDa to 4 MDa. We choose the 64 kDa hemoglobin^{36}, the 706 kDa 20S proteasome from yeast^{2}, and the 4 MDa human ribosome^{37}. We create atomic potential maps using the Matlab code InSilicoTem^{38}, with a thickness of 50 nm and at an electron energy of 300 keV. We use the isolated atom superposition approximation, without solving the PoissonBoltzmann equations for the interaction between the molecule and the ions. We also do not model the amorphousness of the solvent, which was performed in ref. 38 using molecular dynamics simulations, but was seen to have a negligible effect at very low doses. As described in ref. 38, we model the imaginary part of the potential via the inelastic mean free path, creating a minimal transmission contrast between the vitreous ice and the protein. Using these potential maps, we simulate a ptychography experiment by cropping threedimensional slices from the potential at several positions and propagate a coherent incoming wave through the slices using the methods described in ref. 39 in custom code. The final model for the formation of the intensity on the detector is then
for the diffraction pattern and
for the cryoEM image, where \( {\mathcal F} \) and \({ {\mathcal F} }^{1}\) denote the forward and inverse Fourier transform, and CTF the contrast transfer function, respectively.
The detector properties are modeled as in ref. 38, by multiplying the Fourier transform of the exitwave intensity with the square root of the detective quantum efficiency \(\sqrt{{\rm{DQE}}(\vec{q})}\) ^{40}, before applying shot noise by sampling from a Poissonian distribution. We finally convolve the noisy signal with the noise transfer function \({\rm{NTF}}(\vec{q})\) to yield the measured intensity.
A notable difference both in simulation and practice is the fact that for cryoEM, usually no pixel binning is applied to maximize the imaged area and increase throughput. Therefore, also high spatialfrequency regions with low values of DQE and NTF are used for image formation^{41}. For ptychography, on the other hand, the detector can be heavily binned, as long as the realspace patch given by \(\lambda {\rm{\Delta }}z/{d}_{{\rm{pix}}}\equiv {r}_{{\rm{d}}}\cdot {N}_{{\rm{pix}}}\) still encompasses the probe beam on the sample. For typical detectors this condition is fulfilled at bin sizes equivalent to a few percent of the Nyquist frequency. This leads to a nearconstant DQE and a nearunity NTF, such that they can be omitted in the ptychography reconstructions, whereas we still include them in the simulation of the diffraction data. We note, however, that a convolution with a detector transfer function can be modeled with a partially coherent beam if necessary, as demonstrated in refs 42 and 43. We choose the Gatan K2 Summit as the detector for our simulations because it has the highest published DQE and MTF values at low spatial frequencies at 300 keV^{41}. We note that direct detection cameras with frame rates of 1 kHz and above may be more suitable for highthroughput scanning experiments^{44,45,46}, but characteristics for these cameras at 300 keV are either not published or inferior to the K2 Summit. Assuming the K2 Summit for both ptychography and phasecontrast TEM simulations also simplifies a direct comparison between the two methods. The intensity after detection is modeled as^{38}:
where NTF and DQE are properties of the detector^{40, 45} and Poisson(x) samples from a Poisson distribution with mean x.
Singleparticle reconstruction
Figure 2 shows a comparison of lowdose ptychography reconstructions with currently used methods for singleparticle imaging with electrons: defocusbased cryoEM, and Zernike phase contrast cryoEM with a Volta phaseplate. We choose exemplary doses of 5 e ^{−}/Å^{2} as the typical threshold where the highest resolution details are destroyed^{47} and 20 e ^{−}/Å^{2} as a typical dose at which experiments are performed. We have reversed the contrast in the cryoEM images to simplify the visual comparison with the ptychography reconstructions. To quantitatively assess the image quality, we have computed the 2D Fourier Ring Correlation (FRC)^{48} with the ground truth for the both ptychographic reconstruction and simulated cryoEM images of the macromolecules, as shown in Fig. 3. As ground truth for the images we use the electron counts in a noiseless, aberrationfree phaseplate image. Using the 1bit criterion as a resolution threshold^{48}, the achieved resolutions at 5 e ^{−}/Å^{2} and 20 e ^{−}/Å^{2}, respectively, are 12 Å and 8.9 Å for hemoglobin; 10.9 Å and 9.1 Å for 20S proteasome; and 10.3 Å and 5.4 Å for human ribosome. In the case 20S proteasome, these values are identical to the FRC threshold for the phase plate image; for hemoglobin and human ribosome, the phase plate image yields a slightly better resolution of 8.7 Å and 5.1 Å respectively at a dose of 20 e ^{−}/Å^{2} and 10 Å at a dose of 5 e ^{−}/Å^{2}. As the FRC is insensitive to very small and very large values of signaltonoise ratio (SNR), we also show the spatialfrequency resolved SNR in Fig. 3(d)–(f). We define the SNR as
The SNR of the ptychographic reconstruction is significantly lower than the SNR of the phaseplate image for all three particles at spatial frequencies below 0.1 Å^{−1}. It does however scale better to high spatial frequencies, and does not drop below −15 dB for all particles and all resolution up to nearly 0.5 Å^{−1}, where the SNR is close to two orders of magnitude better than the phaseplate image SNR. This helps ptychography perform better when multiple reconstructions are averaged, because a positive singledigit SNR can be reached with fewer particles.
Effect of averaging
In singleparticle cryoEM, a threedimensional reconstruction can be obtained by collecting a large ensemble of 2D images, before orienting and averaging them in three dimensions, such that a satisfactory SNR is achieved. A similar 3D reconstruction from ptychographic data is out of the scope of this paper. A straightforward approach would be to use the reconstructed 2D phase images as an input to existing 3D reconstruction algorithms^{4, 49}, as is done routinely in ptychographic Xray tomography^{50}; however, it is likely that a dedicated algorithm which reconstructs the 3D structure directly from the diffraction patterns will achieve the best results. Coarse orientation alignment for such an algorithm could be done in real space from 2D reconstructions shown here. We leave the evaluation of such approaches to future studies, and for now concentrate on the achievable 2D resolution when averaging a larger ensemble of particles.
To give a rough estimate how the resolution and SNR achieved by our algorithm scales with averaging over multiple datasets, we simulate 30 independent datasets with identically oriented particles, and average the reconstruction results. We found that the resolution corresponding to the diffraction limit as defined by the probeforming aperture is achieved with roughly 40 averaged reconstructions (see the Supplementary Information). Superresolution beyond this point is in principle possible, but due to the low electron counts and the low contrast transfer for higher angles, more images are needed for further improvements. To limit the amount of necessary computation we use an average of 30 images, where a resolution of roughly 3 Å is reached, close to the resolution corresponding to the probeforming aperture of 2.7 Å. Figure 4(a)–(f) show images of the averaged reconstructions of the three samples, at doses of 5 e ^{−}/Å^{2} and 20 e ^{−}/Å^{2} respectively. We also compare the FRC between respectively 30 averaged reconstructions of 60 independently created ptychographic data sets, to give a resolution estimate. We use here the 1/2bit resolution threshold discussed in ref. 48, which gives a slightly more conservative estimate than the 0.143criterion commonly used for the 3D Fourier Shell correlation in averaged reconstructions for 3D cryoEM. With averaging, a resolution of 3.4 Å is achieved for hemoglobin, 3.1 Å for 20S proteasome and 2.9 Å for human ribosome, all at a dose of 20 e ^{−}/Å^{2}. This shows that cryoelectron ptychography can recover atomic resolution information in 2D from only tens of averaged images, facilitated by the favorable scaling of the SNR with spatial frequency, as discussed above. It may therefore be possible to drastically reduce the number of particles that is required for a successful 3D reconstruction at atomic resolution. A rough tomographic estimate^{51} suggests that the ribosome could be reconstructed to 3 Å resolution in 3D with less than 10000 particles.
Probe and dose dependence
It is wellknown that the phase profile of the ptychographic probe can heavily influence the reconstruction quality^{52,53,54,55,56}. Here, we numerically test three different probes, depicted in Fig. 5, and their influence on the reconstruction SNR at low and high doses: (1) a standard defocused probe with defocus aberration of 400 nm, (2) a defocused Fresnel Zone Plate (FZP), and (3) a randomized probe generated by a holographic phase plate and a conventional lens. Figure 5 depicts these probe in real and Fourier space, and typical diffraction patterns at infinite dose and low dose. The FZP was recently suggested as a phase modulator for brightfield STEM^{57}, because its simple phase modulation allows analytical retrieval of linear phase contrast. However, diffractive optics typically have imperfections due to the manufacturing process, which introduce errors and dose inefficiency if the phase modulation is obtained by a simple fitting procedure. Iterative ptychography algorithms allow for the simultaneous retrieval of the probe wave function^{19, 20}, and therefore offer full flexibility in the design of the phase profile. Empirically, probes with a diffuse phase profile result in better reconstructions; therefore, we test as a third probe a random illumination generated by a holographic phase plate and a focusing lens.
Figure 6 shows the SNR of the three proposed probes as a function of spatial frequency for doses of 20 e ^{−}/Å^{2} and 80 e ^{−}/Å^{2}. It can be seen that the simple defocused probe has almost 2 orders of magnitude worse SNR than the FZP and the random probe at the lowest spatial frequencies. At low spatial frequency the FZP achieves the best SNR, while at high spatial frequencies the random probe does slightly better. We have therefore used the random probe for the reconstructions shown in Figs 2 and 4. We give a qualitative explanation of this fact, but emphasize that a theory for optimal measurement design in ptychography and a practically feasible implementation of it is still outstanding, and may drastically improve upon the results presented here. From Fig. 5 column c) it can be seen, that the randomized probe induces the strongest intensity fluctuations (speckle) in the diffraction pattern. These very localized fluctuations vary strongly when scanning the sample, while the diffraction pattern from the defocused probe has slowlyvarying features, which correlate strongly with the real space transmission function and vary only weakly when scanning the sample, which leads to a less diverse dataset and can cause problems in the reconstruction.
Methods
Mathematical framework of Ptychography
We define the twodimensional grid with size \({n}_{1}\times {n}_{2}\in {\mathbb{N}}\times {\mathbb{N}}\) and length scale r > 0 as \({D}_{r}^{{n}_{1}\times {n}_{2}}\,:={(r\alpha ,r\beta )}_{\alpha ,\beta =0}^{{n}_{1},{n}_{2}}\subset {{\mathbb{R}}}^{2}\). The twodimensional complex transmission function of the object is discretized as a n _{1} × n _{2} matrix and denoted as \(T:{D}_{{r}_{{\rm{d}}}}^{{n}_{1}\times {n}_{2}}\to {\mathbb{C}}\), where r _{d} > 0 is the diffractionlimited length scale as introduced above. The object is illuminated by a small beam with known distribution, and discretized as a m _{1} × m _{2} matrix, denoted as \(\psi :{D}_{{r}_{{\rm{d}}}}^{{m}_{1}\times {m}_{2}}\to {\mathbb{C}}\). For simplicity, in this paper we only consider the case n _{1} = n _{2} and m _{1} = m _{2}, i.e. a uniform discretization in both axes. In the experiment, the beam is moved over the sample to positions \({\vec{r}}_{i}\), and illuminates K > 1 subregions to obtain K diffraction images. The intensity measured for position i is then
where the realspace coordinates are discretized in steps of r _{d}, and reciprocalspace coordinates in steps of (m _{{1,2}} r _{d})^{−1}. Mathematically, ptychographic reconstruction can be understood as a special case of the generalized phase retrieval problem: given a phaseless vector of measurements \({\bf{y}}\in {{\mathbb{R}}}_{+}^{m}\) find a complex vector \({\bf{z}}\in {{\mathbb{C}}}^{n}\) such that
where \({\mathscr{A}}:{{\mathbb{C}}}^{n}\to {{\mathbb{C}}}^{m}\) is an arbitrary linear operator. We follow the notations in ref. 52 to write the ptychography problem in this form. First, we vectorize the transmission function as \({{\bf{T}}}^{V}\in {{\mathbb{C}}}^{N}\) with \(N={n}_{1}\cdot {n}_{2}\in {\mathbb{N}}\). We introduce the matrix \({R}_{(i)}\in {{\mathbb{R}}}^{M\times N}\), which extracts an m _{1} × m _{2} sized area centered at position \({\vec{r}}_{i}\) from T ^{V}. With these notations in place, the relation between the noisefree diffraction measurements collected in a ptychography experiment and T ^{V} can be represented compactly as
where P is constructed by cropping K regions from T ^{V} and multiplying by the incoming beam in Q, and applying a 2D discrete Fourier transform F, i.e. P = FQ:
The matrix \({\bf{P}}\in {{\mathbb{C}}}^{KM\times N}\) is sometimes called design matrix, because its entries determine the measurement outcome and reflect the experimental design. In the last decades many algorithms to solve this problem have been devised, only a few of which we will review with regards to lowdose reconstruction in the following section. For the subsequent analysis, we denote the KM row vectors of P as p _{i}.
Bayesian optimization with truncated gradients
The most prominent iterative algorithms to solve the ptychographic phase retrieval problems are the difference map (DM) algorithm^{14}, and the extended ptychographic iterative engine (ePIE)^{20}. The difference map belongs to the family of algorithms which use projections onto nonconvex sets to reach a fixpoint, i.e., the solution lying at the intersection of the two sets. It can be shown that the standard algorithm of alternating projections is equivalent to steepestdescent optimization with a Gaussian likelihood, and is not suited for lowdose reconstructions^{52}, because the for this case Poisson distribution arising from discretized count events differs too strongly from a Gaussian. While this argument does not hold for the more elaborate projection algorithms like DM and relaxed averaged alternating reflections (RAAR)^{58}, they also fail in practice at low doses^{25, 59}, and statistical reconstruction methods have to be used. Thibault and GuizarSicairos^{24} have analyzed maximum likelihood methods in conjunction with a conjugate gradient update rule as a refinement step, after the DM algorithm has converged. They demonstrate improved SNR by two orders of magnitude compared to the DM algorithm alone. They note, however, that starting directly with maximum likelihood optimization often poses convergence problems.
Due to the lack of algorithms with convergence guarantees, the mathematical community has recently picked up the problem, and a host of new algorithms with provable convergence has been developed. While we do not elaborate on them here we point the interested reader to the summary articles^{60, 61} and the article^{62}, which refers to the most recent developments.
Here, we focus on developments which specifically target lowdose applications. Notable in this area is the work by Katkovnik et al.^{63}, which in addition to the maximum likelihood estimate introduces a transformdomain sparsity constraint on the object and optimizes two objective functions in an alternating fashion: one for the maximizing the likelihood, and one for obtaining a sparse representation of the transmission function. However, instead of including the Poissonian likelihood directly, an observation filtering step is performed with a Gaussian likelihood. To obtain a sparse representation of the object, the popular BM3D denoising filter is used^{64}.
During the writing of this paper, Yang et al. suggested using the Wigner Distribution Deconvolution technique for lowdose ptychography^{65}, however no statistical treatment of the measurement process is included so far.
In this work, we formulate ptychographic phase retrieval as a Bayesian inference problem, by rewriting the probability of the transmission function T ^{V} given a set of measurements \({\bf{y}}={({y}_{1},{y}_{2},\ldots ,{y}_{KM})}^{T}\in {{\mathbb{R}}}_{+}^{KM}\) according to Bayes’ rule:
Since the measurements y _{ i } follow the Poissonian distribution
the likelihood function is given by
The prior distribution P(T ^{V}) is usually chosen such that it favors realistic solutions, so that noise is suppressed in the reconstructed image. Here we evaluate two different models. A simple prior, suggested in ref. 42, penalizes large gradients in the image with a Gaussian distribution on the gradient of the transmission function, which is also known as Tikhonov regularization:
with \(\kappa =8\frac{{N}_{pix}^{2}}{{N}_{{\rm{m}}}{\Vert I\Vert }_{1}}\) chosen as in ref. 42. N _{ m } is the total number of valid measurements, N _{ pix } * K in the case when the detector has no hot pixels. This scales the numerical value of the prior to be close to the likelihood, such that the weight μ _{0} can take values between 1 × 10^{−1} and 1 × 10^{−2}. D_{ x } and D_{ y } are the discrete forward difference operators. The second prior we evaluate is based on the work by Katkovnik et al.^{63} and uses sparse modeling to denoise the transmission function:
Here, \({{\bf{T}}}_{{\rm{sparse}}}^{V}\) is built up by applying the BM3D collaborative filtering algorithm^{64, 66}, which we briefly describe in the Supplementary Material. As input for the BM3D algorithm we transform T ^{V} into huesaturationvalue format using domain coloring. The prior P _{sparse}(T) reduces the difference between the denoised version of the current transmission function and the transmission function itself. We do not take into account the marginal likelihood P(y) due to the high dimensionality of the problem. Given the likelihood function P(yT ^{V}) and the prior distribution P(T ^{V}), we can now write the objective function for the maximumaposteriori (MAP) estimate:
The loglikelihood is given as
and the MAP objective functions are
and
for the two prior models, respectively. We calculate the gradients of both expressions:
Since equations (16) and (17) are nonconvex functions, there is no guarantee that standard gradient descent converges to a global minimum. Recently, a nonconvex algorithm for the generalized phase retrieval problem with Poisson noise was presented^{67}, that provably converges to a global minimum with suitable initialization. It introduces a iterationdependent regularization on the gradients of the likelihood to remove terms which have a negative effect on the search direction. Namely, it introduces a truncation criterion
that acts on the gradient of the likelihood and suppresses the influence of measurements that are too incompatible with the reconstruction. The truncation parameter α _{ h } ≥ 5 is described in ref. 67. The regularized likelihood gradient is then
We compute the next step using conjugate gradient descent^{68, 69}, since this lead to much faster convergence compared to the update procedure described in ref. 67.
Initialization
Truncated spectral initialization for ptychography was first proposed by Marchesini et al.^{52}, based on the notion that the highest intensities in the diffraction pattern carry the strongest phase information. They compute the phase of the largest eigenvector of the following hermitian operator:
where ε is chosen such that the largest 20 percent of the intensities are allowed to contribute and F and Q are defined as above. The largest eigenvalue of a sparse hermitian matrix can be efficiently computed either with power iterations^{70}, or with the Arnoldi method^{71}. In ref. 67, truncated spectral initialization with a truncation rule with \({{\bf{1}}}_{{y}_{i} < {\alpha }_{0}^{2}{\lambda }_{0}^{2}}\) is used, with \({\lambda }_{0}=\sqrt{{\sum }_{i=1}^{KM}{y}_{i}}\) and α _{0} a free parameter. It is also important to note that the truncated spectral initialization only produces visually correct initial phase to a dose of roughly 100 e ^{−}/Å^{2}. Figure 7(b) shows an example initialization for a dose of 100 e ^{−}/Å^{2}. For doses below this value, we initialized the transmission function with unity transmission and normaldistributed phase with mean 0.1 and variance of 0.1. With this random initialization we found no problem of convergence for all algorithms tested in this paper.
Reconstruction parameters
All ptychography reconstructions were performed with a probe area overlap of 75% in realspace, where the probe area is defined by all pixels contributing more than 1% of the maximum intensity. This corresponds to a step size of roughly 3 nm, depending on the probe used. At a dose of 20 e ^{−}/Å^{2}, this corresponds to 8540 electrons per diffraction pattern, i.e. an average electron count of 0.52 electrons per detector pixel, and 58 electrons per image pixel for a pixel size of 1.7 Å. For the regularization parameters we performed a grid search evaluating the final normalized root mean square error and found the best values to be μ _{0} = 1 × 10^{−2}, μ _{1} = 8 × 10^{−2}. We choose the biorthogonal spline wavelet transform as the linear transform for BM3D as it achieves the best PSNR for high noise^{72}. For the figures shown in the paper we use the BM3D denoiser. A typical intermediate denoised image for ribosome after 60 iterations can be seen in Fig. 7(c).
Implementation Details
The algorithms presented in this paper were implemented with the Torch scientific computing framework^{73}. The gradient update routines were adapted from the optim package for Torch^{69}. For efficient computing on the graphics processing unit (GPU) with complex numbers, the zcutorch library for CUDA was developed^{74}. Hyperparameter optimization was done with the hypero^{75} package for Torch. For BM3D denoising we use the C++ implementation^{76}. The code was run on an Intel i76700 processor with 32 GB RAM and a NVidia Titan X GPU with 12 GB RAM. The run time for optimization with \({ {\mathcal L} }_{{\rm{MAP}}}\) was 26 s, and for optimization with \({ {\mathcal L} }_{{\rm{BM3D}}{\rm{MAP}}}\) 35 s. This is expected, because the BM3D algorithm used here is not yet implemented on the GPU, and the BM3D denoising is computationally more intensive.
Data Availability
The datasets generated and analysed during the current study are available from the corresponding author on reasonable request.
Conclusion
In this paper we have, via numerical experiments, demonstrated the feasibility of retrieving highresolution electron transmission phase information of biological macromolecules using ptychography and Bayesian optimization. With the methods presented in this paper, it should be possible to achieve a resolution better than 1 nm for true singleparticle imaging of molecular complexes with molecular weights ranging from below 100 kDa to a few MDa, and a resolution around 3 Å with simple averaging of 30 datasets. We have given a detailed explanation of the optimization and initialization procedures used, and have emphasized the importance of choosing an appropriate illumination function. We note that, while the high data redundancy in a ptychographic dataset empirically makes it experimentally very robust, there is much room for improvement in terms of measurement complexity. For the results presented here, the measurement dimension KM is larger than the problem dimension N by a factor of at least 30, while the theoretical limit for successful phase retrieval is at KM = 4N ^{77}. By reducing the number of measurements, the SNR of each individual diffraction measurement could be increased, yielding an improved image SNR in the reconstruction. Therefore, the development of an optimized experimental scheme, including design of the illumination function and scanning scheme is a promising direction of research and may enable significant improvement to the results presented here.
We would like to point out two obstacles that one may have to overcome in the experimental implementation of our method. Firstly, the best results are to be expected when recording zeroloss diffraction patterns with the use of an energy filter. The energy filter may introduce phase distortions into the diffraction patterns, which may need to be accounted for in the reconstruction algorithm. Secondly, although beaminduced movements are expected to be reduced by a large amount due to spotscanning, the remaining movement may cause problems in the reconstruction. Statistically stationary sample movements can be accounted for in the reconstruction algorithm^{25, 78}, but beaminduced motions are likely to be nonstationary, and dedicated algorithms may need to be developed to account for it.
Cryogenic ptychographic imaging of biological samples is also being developed in the Xray sciences^{79}, and our results could equally be implemented there to improve the doseeffectiveness. Finally, the methods presented here may find application in electron phase imaging of radiationsensitive samples under noncryogenic conditions, and the incorporation of Bayesian methods into infocus ptychographic reconstruction procedures^{18, 65}, may provide similar gains in SNR as the ones discussed here, while also keeping the analytical capabilities of traditional scanning TEM imaging.
References
 1.
Bartesaghi, A., Matthies, D., Banerjee, S., Merk, A. & Subramaniam, S. Structure of βgalactosidase at 3.2Å resolution obtained by cryoelectron microscopy. PNAS 111, 11709–11714, doi:10.1073/pnas.1402809111 (2014).
 2.
Groll, M. et al. Structure of 20S proteasome from yeast at 2.4 Å resolution. Nature 386, 463–471, doi:10.1038/386463a0 (1997).
 3.
Khoshouei, M., Radjainia, M., Baumeister, W. & Danev, R. Cryoem structure of haemoglobin at 3.2 a determined with the volta phase plate. Nature Communications, doi:10.1038/ncomms16099 (2017).
 4.
Scheres, S. H. W. RELION: Implementation of a Bayesian approach to cryoEM structure determination. J. Struct. Biol. 180, 519–530, doi:10.1016/j.jsb.2012.09.006 (2012).
 5.
Bai, X.c, Fernandez, I. S., McMullan, G. & Scheres, S. H. Ribosome structures to nearatomic resolution from thirty thousand cryoEM particles. eLife 2, e00461, doi:10.7554/eLife.00461 (2013).
 6.
Zheng, S. Q. et al. MotionCor2: Anisotropic correction of beaminduced motion for improved cryoelectron microscopy. Nat Meth 14, 331–332, doi:10.1038/nmeth.4193 00005 (2017).
 7.
Subramaniam, S., Kühlbrandt, W. & Henderson, R. CryoEM at IUCrJ: A new era. IUCrJ 3, 3–7, doi:10.1107/S2052252515023738 (2016).
 8.
Bai, X.c, McMullan, G. & Scheres, S. H. W. How cryoEM is revolutionizing structural biology. Trends in Biochem. Sci. 40, 49–57, doi:10.1016/j.tibs.2014.10.005 (2015).
 9.
Scheres, S. H. Beaminduced motion correction for submegadalton cryoEM particles. eLife 3, e03665, doi:10.7554/eLife.03665 (2014).
 10.
Glaeser, R. M. How good can cryoEM become? Nat. Meth. 13, 28–32, doi:10.1038/nmeth.3695 (2016).
 11.
Danev, R. & Baumeister, W. CryoEM single particle analysis with the Volta phase plate. eLife 5, e13046, doi:10.7554/eLife.13046 (2016).
 12.
Danev, R., Buijsse, B., Khoshouei, M., Plitzko, J. M. & Baumeister, W. Volta potential phase plate for infocus phase contrast transmission electron microscopy. PNAS 111, 15635–15640, doi:10.1073/pnas.1418377111 (2014).
 13.
Zernike, F. How I Discovered Phase Contrast. Science 121, 345–349, doi:10.1126/science.121.3141.345 (1955).
 14.
Thibault, P. et al. Highresolution scanning xray diffraction microscopy. Science 321, 379–82, doi:10.1126/science.1158573 (2008).
 15.
Maiden, A. M., Sarahan, M. C., Stagg, M. D., Schramm, S. M. & Humphry, M. J. Quantitative electron phase imaging with high sensitivity and an unlimited field of view. Scientific Reports 5, 14690, doi:10.1038/srep14690 (2015).
 16.
Diaz, A. et al. Quantitative xray phase nanotomography. Phys. Rev. B 85, 1–4, doi:10.1103/PhysRevB.85.020104 (2012).
 17.
Hoppe, W. Trace structure analysis, ptychography, phase tomography. Ultramicroscopy 10, 187–198, doi:10.1016/03043991(82)900389 (1982).
 18.
Rodenburg, J. M. The phase problem, microdiffraction and wavelengthlimited resolution — a discussion. Ultramicroscopy 27, 413–422, doi:10.1016/03043991(89)900090 (1989).
 19.
Thibault, P., Dierolf, M., Bunk, O., Menzel, A. & Pfeiffer, F. Probe retrieval in ptychographic coherent diffractive imaging. Ultramicroscopy 109, 338–43, doi:10.1016/j.ultramic.2008.12.011 (2009).
 20.
Maiden, A. M. & Rodenburg, J. M. An improved ptychographical phase retrieval algorithm for diffractive imaging. Ultramicroscopy 109, 1256–62, doi:10.1016/j.ultramic.2009.05.012 (2009).
 21.
GuizarSicairos, M. & Fienup, J. R. Phase retrieval with transverse translation diversity: A nonlinear optimization approach. Opt. Express, OE 16, 7264–7278, doi:10.1364/OE.16.007264 (2008).
 22.
Zhang, F. et al. Translation position determination in ptychographic coherent diffraction imaging. Opt. Express 21, 13592–13606, doi:10.1364/OE.21.013592 (2013).
 23.
Maiden, A. M., Humphry, M. J., Sarahan, M. C., Kraus, B. & Rodenburg, J. M. An annealing algorithm to correct positioning errors in ptychography. Ultramicroscopy 120, 64–72, doi:10.1016/j.ultramic.2012.06.001 (2012).
 24.
Thibault, P. & GuizarSicairos, M. Maximumlikelihood refinement for coherent diffractive imaging. New J. Phys. 14, 063004, doi:10.1088/13672630/14/6/063004 (2012).
 25.
Pelz, P. M. et al. Onthefly scans for Xray ptychography. Appl. Phys. Lett. 105, 251101, doi:10.1063/1.4904943 (2014).
 26.
Clark, J. N., Huang, X., Harder, R. J. & Robinson, I. K. Dynamic Imaging Using Ptychography. Phys. Rev. Lett. 112, 113901, doi:10.1103/PhysRevLett.112.113901 (2014).
 27.
Marchesini, S., Schirotzek, A., Yang, C., Wu, H.t. & Maia, F. Augmented projections for ptychographic imaging. Inverse Problems 29, 115009, doi:10.1088/02665611/29/11/115009 (2013).
 28.
Putkunz, C. T. et al. AtomScale Ptychographic Electron Diffractive Imaging of Boron Nitride Cones. Phys. Rev. Lett. 108, 073901, doi:10.1103/PhysRevLett.108.073901 (2012).
 29.
Yang, H. et al. Simultaneous atomicresolution electron ptychography and Zcontrast imaging of light and heavy elements in complex nanostructures. Nat. Commun. 7, 12532, doi:10.1038/ncomms12532 (2016).
 30.
D’Alfonso, A. J., Allen, L. J., Sawada, H. & Kirkland, A. I. Dosedependent highresolution electron ptychography. J. Appl. Phys. 119, 054302, doi:10.1063/1.4941269 (2016).
 31.
Lubk, A. & Röder, F. Phasespace foundations of electron holography. Phys. Rev. A 92, 033844, doi:10.1103/PhysRevA.92.033844 (2015).
 32.
Humphry, M., Kraus, B. & Hurst, A. Ptychographic electron microscopy using highangle darkfield scattering for subnanometre resolution imaging. Nat. Commun. 3, 730–737, doi:10.1038/ncomms1733 (2012).
 33.
Bullough, P. & Henderson, R. Use of spotscan procedure for recording lowdose micrographs of beamsensitive specimens. Ultramicroscopy 21, 223–230, doi:10.1016/03043991(87)901471 (1987).
 34.
Brink, J., Chiu, W. & Dougherty, M. Computercontrolled spotscan imaging of crotoxin complex crystals with 400 keV electrons at nearatomic resolution. Ultramicroscopy 46, 229–240, doi:10.1016/03043991(92)90017E (1992).
 35.
Downing, K. H. Spotscan imaging in transmission electron microscopy. Science 251, 53–59, doi:10.1126/science.1846047 (1991).
 36.
Fermi, G., Perutz, M. F., Shaanan, B. & Fourme, R. The crystal structure of human deoxyhaemoglobin at 1.74A resolution. J. Mol. Biol. 175, 159–174 (1984).
 37.
Anger, A. M. et al. Structures of the human and Drosophila 80S ribosome. Nature 497, 80–85, doi:10.1038/nature12104 (2013).
 38.
Vulović, M. et al. Image formation modeling in cryoelectron microscopy. J. Struct. Biol. 183, 19–32, doi:10.1016/j.jsb.2013.05.008 (2013).
 39.
Kirkland, E. Advanced Computing in Electron Microscopy (Springer, 2010).
 40.
Meyer, R. R. & Kirkland, A. I. Characterisation of the signal and noise transfer of CCD cameras for electron detection. Microsc. Res. Tech. 49, 269–280, doi:10.1002/(SICI)10970029(20000501)49:3¡269::AIDJEMT5¿3.0.CO;2B (2000).
 41.
McMullan, G., Faruqi, A. R., Clare, D. & Henderson, R. Comparison of optimal performance at 300 keV of three direct electron detectors for use in low dose electron microscopy. Ultramicroscopy 147, 156–163, doi:10.1016/j.ultramic.2014.08.002 (2014).
 42.
Thibault, P. & Menzel, A. Reconstructing state mixtures from diffraction measurements. Nature 494, 68–71, doi:10.1038/nature11806 (2013).
 43.
Enders, B. et al. Ptychography with broadbandwidth radiation. Appl. Phys. Lett. 104, 171104, doi:10.1063/1.4874304 (2014).
 44.
Ryll, H. et al. A pnCCDbased, fast direct single electron imaging camera for TEM and STEM. J. Inst. 11, 4006, doi:10.1088/17480221/11/04/P04006 (2016).
 45.
McMullan, G., Chen, S., Henderson, R. & Faruqi, A. R. Detective quantum efficiency of electron area detectors in electron microscopy. Ultramicroscopy 109, 1126–1143, doi:10.1016/j.ultramic.2009.04.002 (2009).
 46.
Tate, M. W. et al. High Dynamic Range Pixel Array Detector for Scanning Transmission Electron Microscopy. Microsc. Microanal. 22, 237–249, doi:10.1017/S1431927615015664.
 47.
Stark, H., Zemlin, F. & Boettcher, C. Electron radiation damage to protein crystals of bacteriorhodopsin at different temperatures. Ultramicroscopy 63, 75–79, doi:10.1016/03043991(96)000459 (1996).
 48.
Heel, M. V. & Schatz, M. Fourier shell correlation threshold criteria. J. Struct. Biol. 151, 250–262, doi:10.1016/j.jsb.2005.05.009 (2005).
 49.
Tang, G. et al. Eman2: An extensible image processing suite for electron microscopy. J. Struct. Biol. 157, 38–46, doi:10.1016/j.jsb.2006.05.009 (2007).
 50.
Dierolf, M. et al. Ptychographic Xray computed tomography at the nanoscale. Nature 467, 436–439, doi:10.1038/nature09419 (2010).
 51.
Crowther, R. A., DeRosier, D. J. & Klug, A. The reconstruction of a threedimensional structure from projections and its application to electron microscopy. Proceedings of the Royal Society of London. Series A, Mathematical and Physical Sciences 317, 319–340 (1970).
 52.
Marchesini, S., Tu, Y. & Wu, H.t. Alternating Projection, Ptychographic Imaging and Phase Synchronization. arXiv 1402.0550 (2014).
 53.
Li, P. et al. Multiple mode xray ptychography using a lens and a fixed diffuser optic. J. Opt. 18, 054008, doi:10.1088/20408978/18/5/054008 (2016).
 54.
Maiden, A. M., Morrison, G. R., Kaulich, B., Gianoncelli, A. & Rodenburg, J. M. Soft Xray spectromicroscopy using ptychography with randomly phased illumination. Nat. Commun. 4, 1669, doi:10.1038/ncomms2640 (2013).
 55.
GuizarSicairos, M. et al. Role of the illumination spatialfrequency spectrum for ptychography. Phys. Rev. B 86, 100103, doi:10.1103/PhysRevB.86.100103 (2012).
 56.
Li, P., Edo, T. B. & Rodenburg, J. M. Ptychographic inversion via Wigner distribution deconvolution: Noise suppression and probe design. Ultramicroscopy 147, 106–113, doi:10.1016/j.ultramic.2014.07.004 (2014).
 57.
Ophus, C. et al. Efficient linear phase contrast in scanning transmission electron microscopy with matched illumination and detector interferometry. Nat Commun 7, 10719, doi:10.1038/ncomms10719 (2016).
 58.
Luke, D. R. Relaxed averaged alternating reflections for diffraction imaging. Inverse Problems 21, 37, doi:10.1088/02665611/21/1/004 (2005).
 59.
Godard, P., Allain, M., Chamard, V. & Rodenburg, J. Noise models for low counting rate coherent diffraction imaging. Opt. Express 20, 25914–34 (2012).
 60.
Jaganathan, K., Eldar, Y. C. & Hassibi, B. Phase Retrieval: An Overview of Recent Developments. arXiv 1510.07713 (2015).
 61.
Shechtman, Y. et al. Phase Retrieval with Application to Optical Imaging: A contemporary overview. IEEE Signal Processing Magazine 32, 87–109, doi:10.1109/MSP.2014.2352673 (2015).
 62.
Sun, J., Qu, Q. & Wright, J. A geometric analysis of phase retrieval. In 2016 IEEE International Symposium on Information Theory (ISIT), 2379–2383, doi:10.1109/ISIT.2016.7541725 (2016).
 63.
Katkovnik, V. & Astola, J. Sparse ptychographical coherent diffractive imaging from noisy measurements. J. Opt. Soc. Am. A. Opt. Image Sci. Vis. 30, 367–79, doi:10.1364/JOSAA.30.000367 (2013).
 64.
Danielyan, A., Katkovnik, V. & Egiazarian, K. BM3D Frames and Variational Image Deblurring. IEEE Transactions on Image Processing 21, 1715–1728, doi:10.1109/TIP.2011.2176954 (2012).
 65.
Yang, H., Ercius, P., Nellist, P. D. & Ophusa, C. Enhanced phase contrast transfer using ptychography combined with a prespecimen phase plate in a scanning transmission electron microscope. Ultramicroscopy., doi:10.1016/j.ultramic.2016.09.002.
 66.
Chang, H. & Marchesini, S. A general framework for denoising phaseless diffraction measurements. arXiv 1611.01417 (2016).
 67.
Chen, Y. & Candes, E. J. Solving Random Quadratic Systems of Equations Is Nearly as Easy as Solving Linear Systems. arXiv 1505.05114 (2015).
 68.
Shewchuk, J. R. An Introduction to the Conjugate Gradient Method Without the Agonizing Pain. Tech. Rep., Carnegie Mellon University, Pittsburgh, PA, USA (1994).
 69.
torch.optim  numerical optimization library for torch. URL https://github.com/torch/optim/.
 70.
Mises, R. V. & PollaczekGeiringer, H. Praktische Verfahren der Gleichungsauflösung. Z. angew. Math. Mech. 9, 152–164, doi:10.1002/zamm.19290090206 (1929).
 71.
Lehoucq, R., Sorensen, D. & Yang, C. ARPACK Users’ Guide: Solution of LargeScale Eigenvalue Problems with Implicitly Restarted Arnoldi Methods. SIAM ebooks (Society for Industrial and Applied Mathematics (SIAM, 3600 Market Street, Floor 6, Philadelphia, PA 19104), 1998).
 72.
Lebrun, M. A. Analysis and Implementation of the BM3D Image Denoising Method. Image Processing On Line 2, 175–213, doi:10.5201/ipol.2012.lbm3d (2012).
 73.
Collobert, R., Kavukcuoglu, K. & Farabet, C. Torch7: A Matlablike Environment for Machine Learning. In BigLearn, NIPS Workshop (2011).
 74.
Zcutorch  Complex number support for cutorch. https://github.com/PhilippPelz/zcutorch.
 75.
Leonard, N. Hypero  Hyperparameter optimization for torch. https://github.com/ElementResearch/hypero.
 76.
Lebrun, M. bm3d  C++ implementation of BM3D denoising. https://github.com/gfacciol/bm3d.
 77.
Balan, R., Casazza, P. & Edidin, D. On signal reconstruction without phase. Appl. Comput. Harmon. Anal. 20, 345–356, doi:10.1016/j.acha.2005.07.001 (2006).
 78.
Clark, J. N. et al. Dynamic sample imaging in coherent diffractive imaging. Opt. Lett. 36, 1954–6, doi:10.1364/OL.36.001954 (2011).
 79.
Diaz, A. et al. Threedimensional mass density mapping of cellular ultrastructure by ptychographic Xray nanotomography. J. Struct. Biol. 192, 461–469, doi:10.1016/j.jsb.2015.10.008 (2015).
Acknowledgements
This work was funded by the Max Planck Society through institutional support. P.M.P. acknowledges support from the International Max Planck Research School for Ultrafast Imaging & Structural Dynamics.
Author information
Affiliations
Contributions
P.M.P., R.J.D.M. and R.B. conceived the idea. P.M.P. and W.Q. wrote the simulation code. P.M.P. wrote the ptychographic reconstruction code. P.M.P. wrote the paper with input and discussions from R.B., G.K., and R.J.D.M.
Corresponding author
Ethics declarations
Competing Interests
The authors declare that they have no competing interests.
Additional information
Publisher's note: Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Electronic supplementary material
Rights and permissions
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.
About this article
Cite this article
Pelz, P.M., Qiu, W.X., Bücker, R. et al. Lowdose cryo electron ptychography via nonconvex Bayesian optimization. Sci Rep 7, 9883 (2017). https://doi.org/10.1038/s4159801707488y
Received:
Accepted:
Published:
Further reading

A data reduction and compression description for high throughput timeresolved electron microscopy
Nature Communications (2021)

Mixedstate electron ptychography enables subangstrom resolution imaging with picometer precision at low dose
Nature Communications (2020)

Lowdose phase retrieval of biological specimens using cryoelectron ptychography
Nature Communications (2020)

Serial protein crystallography in an electron microscope
Nature Communications (2020)

Fast approximate STEM image simulations from a machine learning model
Advanced Structural and Chemical Imaging (2019)
Comments
By submitting a comment you agree to abide by our Terms and Community Guidelines. If you find something abusive or that does not comply with our terms or guidelines please flag it as inappropriate.