8. Introduction to Cryo-EM#

Cryo-electron microscopy (cryo-EM) is a family of electron microscopy techniques applied to biological specimens preserved by rapid cooling in vitrified — amorphous, non-crystalline — ice at cryogenic temperature. It has become the dominant method for determining the three-dimensional structures of proteins and macromolecular complexes, earning the 2017 Nobel Prize in Chemistry for Jacques Dubochet, Joachim Frank, and Richard Henderson. Unlike X-ray crystallography, cryo-EM does not require crystals and is not limited by molecular size. Unlike NMR, it is sensitive to large and heterogeneous complexes. The key enabling principle is that the electron beam can pass through a thin hydrated specimen that has been vitrified rather than dried or chemically fixed, preserving structures in a state close to their solution conformation.

Two main families of methods share this platform. Single-particle analysis (SPA) averages tens of thousands to millions of images of identical, randomly oriented copies of a molecule to reconstruct a single consensus 3D density map at near-atomic resolution. Electron tomography (cryo-ET) records a series of images at different tilt angles and back-projects them into a 3D volume, giving a unique view of structurally heterogeneous objects such as intact cells, organelles, or in-situ protein assemblies. Both depend on the same physical and mathematical foundations: the physics of electron–specimen interaction, the optics of the transmission electron microscope, and the mathematics of image formation and signal recovery from noise.

This chapter introduces these shared foundations, building on the TEM instrument and contrast transfer function described in Section 7. The following chapters develop the mathematics of Fourier analysis and image formation (Section 9), the SPA workflow (Section 10), and electron tomography (Section 11) in detail.

8.1. The challenge of imaging biological specimens#

Biological macromolecules are built from light elements — carbon, nitrogen, oxygen, and hydrogen — with low nuclear charge. They therefore scatter electrons weakly compared with heavy-metal stains or inorganic materials. Two physical properties combine to make this an unusually difficult imaging problem.

Low contrast. Unstained biological material produces almost no amplitude contrast: it scatters a tiny fraction of electrons (a few percent for a thin protein) and is nearly transparent to the beam. The dominant image-forming mechanism is instead phase contrast — the unscattered beam and the slightly deflected, phase-shifted scattered beam interfere after the objective lens to produce contrast. This mechanism is efficient only for specific spatial frequencies selected by the defocus of the objective lens (the contrast transfer function, or CTF; see Section 7.4).

Radiation damage. Energetic electrons (100–300 keV) deposit energy in the specimen through inelastic scattering, breaking chemical bonds, displacing atoms, and degrading the structure. The total electron dose a biological specimen can sustain before structural information is lost is approximately 40–100 electrons per square ångström (e⁻/Ų). This absolute physical limit means that the number of electrons available to form the image is fundamentally constrained, and there is no instrumental solution that allows “exposing longer” to improve signal.

These two properties together — weak phase contrast and strict dose limits — mean that cryo-EM images are inherently noisy. The strategy for dealing with this noise depends on the imaging modality: SPA averages over many identical copies, exploiting the statistical law that averaging \(N\) images improves the signal-to-noise ratio (SNR) by \(\sqrt{N}\); cryo-ET accepts the noise in individual tilt images and instead uses algorithmic reconstruction to isolate structural information.

8.2. Signal, noise, and the need for averaging#

Electrons are discrete quanta. When \(n\) electrons impinge on a pixel of the detector, the measured count follows a Poisson distribution with variance \(n\). The irreducible shot noise limits the SNR of a single measurement:

\[ \mathrm{SNR}_\mathrm{single} = \frac{\bar{n}}{\sqrt{\bar{n}}} = \sqrt{\bar{n}} \]

In cryo-EM, a typical dose of 50 e⁻/Ų spread over a micrograph with 1 Å pixels delivers about 50 electrons per pixel. Shot noise then contributes \(\mathrm{SNR} \approx \sqrt{50} \approx 7\). But the signal of interest — the contrast modulation produced by a single ~100 kDa protein embedded in ice — contributes only a small fraction of that total, because most of the electrons pass through the surrounding vitreous ice rather than through the protein. The SNR per particle image is typically 0.01–0.1, so the protein is literally invisible in a single raw image.

The solution is coherent averaging. If \(N\) images of identical objects are aligned and averaged, the signal (which is the same in each image) accumulates proportionally to \(N\), while the independent random noise accumulates only as \(\sqrt{N}\):

(8.1)#\[ \mathrm{SNR}(N) = \mathrm{SNR}(1) \times \sqrt{N} \]

To improve SNR by a factor of 10, one needs \(N = 100\) images; a factor of 100 requires \(N = 10^4\). Modern SPA datasets routinely contain \(10^5\)\(10^6\) particle images, achieving SNR improvements of 300–1000× over a single exposure. The interactive below demonstrates how shot noise degrades images at low dose and how averaging recovers the underlying signal.

Interactive element

Click Live Code to activate, expand the Show code toggle, and click ▶ to run. Use the dose slider to set the number of electrons per Ų per image, and observe how the image quality improves as more images are averaged together.

Hide code cell source

import io
import numpy as np
import matplotlib
matplotlib.use('agg')
import matplotlib.pyplot as plt
from ipywidgets import IntSlider, Output, VBox, HBox, Label
from IPython.display import display, Image, clear_output

_rng = np.random.default_rng(7)

def _make_molecule(n=80):
    x = np.linspace(-1, 1, n)
    X, Y = np.meshgrid(x, x)
    p  = 0.9 * np.exp(-(X**2 + Y**2)/(2*0.20**2))
    p += 0.6 * np.exp(-((X-0.35)**2 + (Y+0.05)**2)/(2*0.14**2))
    p += 0.6 * np.exp(-((X+0.32)**2 + (Y-0.10)**2)/(2*0.14**2))
    p += 0.4 * np.exp(-((X+0.05)**2 + (Y+0.40)**2)/(2*0.10**2))
    p += 0.3 * np.exp(-((X-0.15)**2 + (Y-0.42)**2)/(2*0.08**2))
    return p / p.max()

_mol = _make_molecule()

def draw_dose(dose=20):
    epp = dose * 1.3**2
    signal = (_mol * epp).clip(0)
    fig, axes = plt.subplots(1, 4, figsize=(13, 3.5))
    for ax, n_avg in zip(axes, [1, 4, 16, 64]):
        imgs = _rng.poisson(signal, size=(n_avg, 80, 80)).astype(float)
        avg = imgs.mean(0); avg -= avg.mean()
        std = avg.std(); vr = max(std * 3, 0.01)
        ax.imshow(avg, cmap='gray', vmin=-vr, vmax=vr)
        ax.set_title(f"N = {n_avg}\n(SNR ∝ {np.sqrt(n_avg):.1f}×)", fontsize=9)
        ax.axis('off')
    fig.suptitle(f"Dose = {dose} e⁻/Ų per image", fontsize=10)
    fig.tight_layout()
    buf = io.BytesIO(); fig.savefig(buf, format='png', dpi=96); buf.seek(0)
    display(Image(data=buf.read())); plt.close(fig)

style = {"description_width": "130px"}
dose_sl = IntSlider(value=20, min=2, max=100, step=2,
                    description="Dose (e⁻/Ų)", style=style,
                    layout={"width": "380px"})

out = Output()
def update_dose(_=None):
    with out:
        clear_output(wait=True)
        draw_dose(dose_sl.value)
dose_sl.observe(update_dose, names='value')

display(VBox([dose_sl, out]))
update_dose()

8.3. Phase contrast and the Contrast Transfer Function#

As noted above, biological specimens are almost pure phase objects. To understand why this matters, consider the exit wave \(\psi(\mathbf{r})\) after the specimen. For a thin phase object with projected potential \(V(\mathbf{r})\):

\[ \psi(\mathbf{r}) = \exp\!\bigl[i\sigma V(\mathbf{r})\bigr] \approx 1 + i\sigma V(\mathbf{r}) \]

where \(\sigma\) is the electron interaction parameter (of order \(10^{-3}\) V⁻¹ nm for 300 keV electrons) and the approximation holds for weak scatterers (the weak phase object approximation, WPO). The key point is that the second term, \(i\sigma V\), is purely imaginary — it is 90° out of phase with the unscattered beam (the “1” term). The intensity at the exit plane is:

\[ |\psi|^2 = |1 + i\sigma V|^2 \approx 1 + 2\,\mathrm{Im}(\sigma V) + |\sigma V|^2 \approx 1 \]

A pure phase object produces no detectable image contrast because the detector responds only to intensity \(|\psi|^2\), not to phase. This is the fundamental problem: the protein structure is encoded in the phase of the exit wave, but the detector is blind to phase.

Phase contrast is recovered by introducing a phase shift between the unscattered beam and the scattered beam, rotating their relative phase from 90° toward 0° or 180°. When the two beams are in phase (or anti-phase), their interference generates detectable amplitude (intensity) contrast. In a light microscope this is done with a phase plate; in a TEM it is achieved by defocusing the objective lens. Defocus introduces a spatially varying phase shift that depends on the spatial frequency \(s\) of the scattered wave — precisely the contrast transfer function (CTF). The CTF converts phase-object information into amplitude contrast, but does so unevenly: it is a sinusoidal function of \(s\) that passes some frequencies with high contrast, inverts others, and produces zeros (no contrast) at specific spatial frequencies. Understanding and correcting the CTF is a central step in any cryo-EM data processing pipeline (see Section 7.4).

8.4. Fourier analysis in cryo-EM imaging#

The natural language of image formation in the TEM is the Fourier transform. The Fourier transform \(\hat{f}(\mathbf{k})\) of an image \(f(\mathbf{r})\) decomposes it into spatial frequency components, where the spatial frequency vector \(\mathbf{k}\) has magnitude \(s = |\mathbf{k}|\) (units: Å⁻¹) and the corresponding real-space spacing is \(d = 1/s\) (units: Å). Low spatial frequencies encode slowly varying, coarse features; high spatial frequencies encode fine structural detail:

\[ \hat{f}(\mathbf{k}) = \int f(\mathbf{r})\, e^{-2\pi i \mathbf{k} \cdot \mathbf{r}}\, d\mathbf{r} \]

The resolution of a cryo-EM reconstruction is the highest spatial frequency \(s_\mathrm{max}\) at which signal can be reliably measured — equivalently, the smallest feature spacing \(d_\mathrm{min} = 1/s_\mathrm{max}\). Removing spatial frequencies above a cutoff (low-pass filtering) blurs the image and removes fine detail, while the Fourier power spectrum \(|\hat{f}(\mathbf{k})|^2\) shows which frequencies carry the most power.

Three results from Fourier theory are central to cryo-EM:

Convolution theorem. The image recorded by the TEM is the ideal object convolved with the point spread function (PSF) of the microscope. In Fourier space this becomes a pointwise multiplication: \(\hat{I}(\mathbf{k}) = \hat{O}(\mathbf{k}) \cdot \mathrm{CTF}(\mathbf{k})\). CTF correction is therefore a division (or Wiener filter) in Fourier space.

Projection slice theorem. The 2D Fourier transform of a projection of a 3D object is a central 2D slice through the 3D Fourier transform of that object. Collecting projections from many angles fills in 3D Fourier space slice by slice; inverse Fourier transformation then recovers the 3D density. This is the foundation of both SPA (random angle coverage) and tomographic reconstruction (controlled tilt series), discussed in Section 10 and Section 11.

Parseval’s theorem. The total signal power is the same whether measured in real space or Fourier space: \(\sum_\mathbf{r}|f|^2 = \sum_\mathbf{k}|\hat{f}|^2\). This allows SNR and resolution to be assessed equivalently in either domain — the basis of Fourier Shell Correlation (FSC) used for resolution estimation.

The interactive below lets you explore how spatial frequency content determines image detail.

Interactive element

Click Live Code to activate, expand the Show code toggle, and click ▶ to run. Move the resolution cutoff slider and observe how removing high-frequency components blurs the image, and where the signal is concentrated in the Fourier power spectrum.

Hide code cell source

import io
import numpy as np
import matplotlib
matplotlib.use('agg')
import matplotlib.pyplot as plt
from ipywidgets import FloatSlider, Output, VBox
from IPython.display import display, Image, clear_output

def _make_mol(n=80):
    x = np.linspace(-1, 1, n)
    X, Y = np.meshgrid(x, x)
    p  = 0.9 * np.exp(-(X**2 + Y**2)/(2*0.20**2))
    p += 0.6 * np.exp(-((X-0.35)**2 + (Y+0.05)**2)/(2*0.14**2))
    p += 0.6 * np.exp(-((X+0.32)**2 + (Y-0.10)**2)/(2*0.14**2))
    p += 0.4 * np.exp(-((X+0.05)**2 + (Y+0.40)**2)/(2*0.10**2))
    p += 0.3 * np.exp(-((X-0.15)**2 + (Y-0.42)**2)/(2*0.08**2))
    return p / p.max()

_mol2 = _make_mol()
_n2 = _mol2.shape[0]
_pixel_A = 1.3
_freqs = np.fft.fftshift(np.fft.fftfreq(_n2, d=_pixel_A))
_Fx, _Fy = np.meshgrid(_freqs, _freqs)
_R2 = np.sqrt(_Fx**2 + _Fy**2)
_ft2 = np.fft.fftshift(np.fft.fft2(_mol2))
_power2 = np.abs(_ft2)**2

def draw_fourier(res_A=5.0):
    cutoff = 1.0 / res_A
    mask = _R2 <= cutoff
    filtered = np.real(np.fft.ifft2(np.fft.ifftshift(_ft2 * mask)))

    fig, axes = plt.subplots(1, 3, figsize=(12, 4))

    axes[0].imshow(_mol2, cmap='gray')
    axes[0].set_title("Original image", fontsize=10); axes[0].axis('off')

    axes[1].imshow(np.log1p(_power2), cmap='viridis')
    theta = np.linspace(0, 2*np.pi, 300)
    r_px = cutoff / _freqs.max() * _n2/2
    axes[1].plot(_n2/2 + r_px*np.cos(theta), _n2/2 + r_px*np.sin(theta),
                 'r-', linewidth=1.5, label=f'{res_A:.0f} Å cutoff')
    axes[1].set_title("Fourier power spectrum", fontsize=10)
    axes[1].legend(fontsize=8, loc='upper right'); axes[1].axis('off')

    axes[2].imshow(filtered, cmap='gray')
    axes[2].set_title(f"Low-pass filtered\n(resolution ≥ {res_A:.0f} Å)", fontsize=10)
    axes[2].axis('off')

    fig.tight_layout()
    buf = io.BytesIO(); fig.savefig(buf, format='png', dpi=96); buf.seek(0)
    display(Image(data=buf.read())); plt.close(fig)

res_sl = FloatSlider(value=5.0, min=2.0, max=25.0, step=0.5,
                     description="Cutoff (Å)",
                     style={"description_width": "100px"},
                     layout={"width": "380px"})

out2 = Output()
def update_fourier(_=None):
    with out2:
        clear_output(wait=True)
        draw_fourier(res_sl.value)
res_sl.observe(update_fourier, names='value')

display(VBox([res_sl, out2]))
update_fourier()

8.5. Cross-correlation and image alignment#

A recurring problem in cryo-EM is alignment: given a reference image \(A\) and a noisy observed image \(X_i\), find the rotation (and translation) that best superimposes the two. This is solved by the cross-correlation (CC) function.

For two images \(f\) and \(g\), the cross-correlation at displacement \(\mathbf{t}\) is:

\[ \mathrm{CC}(\mathbf{t}) = \sum_{\mathbf{r}} f(\mathbf{r})\, g(\mathbf{r} + \mathbf{t}) \]

The displacement \(\mathbf{t}^*\) that maximises CC is the best estimate of the relative shift between the images. By the convolution theorem, this can be computed efficiently in Fourier space:

\[ \mathrm{CC}(\mathbf{t}) = \mathcal{F}^{-1}\!\left[\hat{f}(\mathbf{k})\,\hat{g}^*(\mathbf{k})\right] \]

making alignment of even large images computationally tractable.

For rotation alignment, the reference \(A\) is rotated through a set of candidate angles \(\{\theta_j\}\), and the CC score between the rotated reference \(R^{\theta_j} A\) and the observed image \(X_i\) is computed for each angle. The angle maximising CC is the best-fit rotation. This is equivalent to minimising the squared difference (SQD):

\[ \mathrm{SQD}(\theta) = \sum_{\mathbf{r}} \bigl[X_i(\mathbf{r}) - R^\theta A(\mathbf{r})\bigr]^2 = \mathrm{const} - 2\,\mathrm{CC}(\theta) \]

where the constant terms (sum of squared pixel values of \(X_i\) and of \(A\)) do not depend on \(\theta\). Minimising SQD and maximising CC are therefore equivalent.

Maximum-likelihood (ML) alignment generalises this further. At the low SNR typical of cryo-EM, the CC peak is broad and noisy: many candidate angles give nearly equal CC scores. Hard assignment of the single best-fit angle ignores this uncertainty and can introduce systematic errors. ML instead treats each orientation as a probability:

\[ P(\theta_j \mid X_i, A) \propto \exp\!\left[\frac{\mathrm{CC}(\theta_j)}{\sigma^2}\right] \]

where \(\sigma^2\) is the noise variance. The reconstruction is then a weighted average over all angles, with each angle contributing in proportion to its likelihood. At high noise, the distribution is broad (many angles contribute equally); at low noise, it collapses to the single best-fit angle, recovering the hard-assignment limit. This is the statistical foundation of RELION, cryoSPARC, and other modern SPA software (see Section 10.6).

The interactive below illustrates how the CC score varies with rotation angle and how noise broadens the peak, making alignment increasingly uncertain.

Interactive element

Click Live Code to activate, expand the Show code toggle, and click ▶ to run. Adjust the true rotation angle and the noise level. Observe how a sharp CC peak at low noise becomes broad and ambiguous at high noise — motivating ML alignment over hard assignment.

Hide code cell source

import io
import numpy as np
import matplotlib
matplotlib.use('agg')
import matplotlib.pyplot as plt
from scipy.ndimage import rotate as nd_rotate
from ipywidgets import FloatSlider, IntSlider, VBox, HBox, Output
from IPython.display import display, Image, clear_output

def _make_ref(n=64):
    x = np.linspace(-1, 1, n)
    X, Y = np.meshgrid(x, x)
    p  = 0.9 * np.exp(-(X**2 + Y**2)/(2*0.25**2))
    p += 0.7 * np.exp(-((X-0.40)**2 + Y**2)/(2*0.15**2))
    p += 0.5 * np.exp(-(X**2 + (Y-0.42)**2)/(2*0.10**2))
    p -= 0.3 * np.exp(-((X+0.30)**2 + (Y+0.30)**2)/(2*0.10**2))
    return p

_ref2 = _make_ref()
_angles_cc = np.arange(0, 360, 2)
_rots_cc = [nd_rotate(_ref2, a, reshape=False, order=1) for a in _angles_cc]
_rng_cc = np.random.default_rng(42)

def _ncc(a, b):
    a0 = a - a.mean(); b0 = b - b.mean()
    d = np.sqrt((a0**2).sum() * (b0**2).sum())
    return (a0 * b0).sum() / d if d > 1e-10 else 0.0

def draw_cc(true_angle=75, noise_sigma=0.5):
    noisy = nd_rotate(_ref2, true_angle, reshape=False, order=1) + \
            _rng_cc.normal(0, noise_sigma, _ref2.shape)
    scores = [_ncc(r, noisy) for r in _rots_cc]
    best = _angles_cc[np.argmax(scores)]

    fig, axes = plt.subplots(1, 3, figsize=(13, 4))

    axes[0].imshow(_ref2, cmap='gray')
    axes[0].set_title("Reference image", fontsize=10); axes[0].axis('off')

    vr = max(abs(noisy.max()), abs(noisy.min()))
    axes[1].imshow(noisy, cmap='gray', vmin=-vr, vmax=vr)
    axes[1].set_title(f"Observed image (rotated {true_angle}°, σ={noise_sigma:.1f})", fontsize=10)
    axes[1].axis('off')

    axes[2].plot(_angles_cc, scores, color='steelblue', linewidth=1.5)
    axes[2].axvline(true_angle % 360, color='green', linestyle='--',
                    linewidth=1.8, label=f'True angle: {true_angle}°')
    axes[2].axvline(best, color='red', linestyle=':',
                    linewidth=1.8, label=f'CC peak: {best}°')
    axes[2].set_xlabel("Candidate rotation (°)"); axes[2].set_ylabel("Normalized CC")
    axes[2].set_title("Cross-correlation vs. rotation angle", fontsize=10)
    axes[2].legend(fontsize=9); axes[2].set_xlim(0, 358)

    fig.tight_layout()
    buf = io.BytesIO(); fig.savefig(buf, format='png', dpi=96); buf.seek(0)
    display(Image(data=buf.read())); plt.close(fig)

style = {"description_width": "130px"}
ang_sl   = IntSlider(value=75, min=0, max=178, step=2,
                     description="True angle (°)", style=style, layout={"width":"360px"})
noise_sl = FloatSlider(value=0.5, min=0.1, max=3.0, step=0.1,
                       description="Noise σ", style=style, layout={"width":"360px"})

out3 = Output()
def update_cc(_=None):
    with out3:
        clear_output(wait=True)
        draw_cc(ang_sl.value, noise_sl.value)
for sl in [ang_sl, noise_sl]:
    sl.observe(update_cc, names='value')

display(VBox([ang_sl, noise_sl, out3]))
update_cc()