Broad-Allen Developing Mouse Brain E15.5: Align¶

This tutorial uses celldega>=0.26.1. It shows how Celldega aligns serial sections into a 3D stack in two simple steps, using the Leiden clusters as landmarks:

1. Cluster centroids: one landmark per (slice, cluster) at the cluster's centroid, then a rigid Procrustes fit chained out from the middle slice.
2. Neighborhoods: one alpha shape per (slice, cluster); each slice's rigid pose is refined to maximize the overlap of matching shapes in nearby slices.

The aligned cells are then written as CellCloud DegaFiles, which the companion notebook Broad-Allen_Dev-Mouse-Brain_E15.5_Viz visualizes.

The experiment¶

This dataset is being generated in collaboration with the Allen Institute. It is whole-brain spatial transcriptomics of a mouse embryo at E15.5. The brain was cut into 44 coronal sections, ordered from the front of the brain (slice 1, anterior) to the back (slice 44, posterior), with a section taken every 100 µm, and profiled with the 10x Xenium platform and a 480-gene panel. Cells were segmented from the images, which gives each cell a position (in µm) and a gene-count vector. After quality control the dataset holds about 3.8 million cells.

All cells from all slices were clustered together with Leiden using rapids-singlecell on Manifold.ai, so a cluster means the same thing in every slice. The clustered AnnData (raw counts, cell positions, slice numbers and Leiden labels) is hosted on Hugging Face: Allen-Broad_Dev-Mouse-Brain/E15.5/xenium_E15.5_clustered.h5ad (about 1.9 GB).

In [1]:

import itertools
from pathlib import Path
import shutil
from urllib.request import Request, urlopen
import warnings

import anndata as ad
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import shapely

import celldega as dega


# anndata notes each per-slice subset it re-indexes; harmless here
warnings.filterwarnings("ignore", message="Transforming to str index")

dega.__version__

import itertools from pathlib import Path import shutil from urllib.request import Request, urlopen import warnings import anndata as ad import matplotlib.pyplot as plt import numpy as np import pandas as pd import shapely import celldega as dega # anndata notes each per-slice subset it re-indexes; harmless here warnings.filterwarnings("ignore", message="Transforming to str index") dega.__version__

Out[1]:

'0.26.1'

Data¶

The clustered AnnData is downloaded from Hugging Face, unless it is already on disk.

In [2]:

HF_RESOLVE_BASE = (
    "https://huggingface.co/datasets/broadinstitute/Celldega_Supporting_Data/resolve/main"
)
H5AD_REMOTE = "Allen-Broad_Dev-Mouse-Brain/E15.5/xenium_E15.5_clustered.h5ad"

DATA_DIR = Path("data/Allen_Dev_Mouse_Brain/E15.5")
H5AD_PATH = DATA_DIR / "xenium_E15.5_clustered.h5ad"


def download_if_missing(remote_filename, destination):
    """Download one public Hub file, with no HF client dependency."""
    if destination.is_file():
        return destination

    destination.parent.mkdir(parents=True, exist_ok=True)
    partial = destination.with_name(f"{destination.name}.part")
    url = f"{HF_RESOLVE_BASE}/{remote_filename}?download=true"
    request = Request(url, headers={"User-Agent": "celldega-example-notebook"})
    print(f"Downloading {url}")
    try:
        with urlopen(request) as response, partial.open("wb") as output:
            shutil.copyfileobj(response, output, length=8 * 1024 * 1024)
        partial.replace(destination)
    except Exception:
        partial.unlink(missing_ok=True)
        raise
    return destination


download_if_missing(H5AD_REMOTE, H5AD_PATH)

HF_RESOLVE_BASE = ( "https://huggingface.co/datasets/broadinstitute/Celldega_Supporting_Data/resolve/main" ) H5AD_REMOTE = "Allen-Broad_Dev-Mouse-Brain/E15.5/xenium_E15.5_clustered.h5ad" DATA_DIR = Path("data/Allen_Dev_Mouse_Brain/E15.5") H5AD_PATH = DATA_DIR / "xenium_E15.5_clustered.h5ad" def download_if_missing(remote_filename, destination): """Download one public Hub file, with no HF client dependency.""" if destination.is_file(): return destination destination.parent.mkdir(parents=True, exist_ok=True) partial = destination.with_name(f"{destination.name}.part") url = f"{HF_RESOLVE_BASE}/{remote_filename}?download=true" request = Request(url, headers={"User-Agent": "celldega-example-notebook"}) print(f"Downloading {url}") try: with urlopen(request) as response, partial.open("wb") as output: shutil.copyfileobj(response, output, length=8 * 1024 * 1024) partial.replace(destination) except Exception: partial.unlink(missing_ok=True) raise return destination download_if_missing(H5AD_REMOTE, H5AD_PATH)

Out[2]:

PosixPath('data/Allen_Dev_Mouse_Brain/E15.5/xenium_E15.5_clustered.h5ad')

In [3]:

CLUSTER_KEY = "leiden_1.0"  # Leiden resolution used throughout
Z_SPACING_UM = 100.0  # distance between consecutive slices in the 3D stack
ALPHA_UM = 100.0  # alpha-shape inv_alpha (circumradius threshold), um
DROP_SLICES = [44]  # slice(s) left out of the 3D stack

# Stage 1, centroid Procrustes
PROCRUSTES_KW = {"alignment_window": 10, "weight_by_adjacent_counts": False}
# Stage 2, neighborhood overlap
NBHD_KW = {
    "alignment_window": 4,
    "n_sweeps": 4,
    "rotation_range": 40.0,
    "translation_range": 0.20,
    "n_rotation_grid": 15,
    "distance_weight": "uniform",
}

CLUSTER_KEY = "leiden_1.0" # Leiden resolution used throughout Z_SPACING_UM = 100.0 # distance between consecutive slices in the 3D stack ALPHA_UM = 100.0 # alpha-shape inv_alpha (circumradius threshold), um DROP_SLICES = [44] # slice(s) left out of the 3D stack # Stage 1, centroid Procrustes PROCRUSTES_KW = {"alignment_window": 10, "weight_by_adjacent_counts": False} # Stage 2, neighborhood overlap NBHD_KW = { "alignment_window": 4, "n_sweeps": 4, "rotation_range": 40.0, "translation_range": 0.20, "n_rotation_grid": 15, "distance_weight": "uniform", }

1. Load the cells¶

Only the cell metadata, positions and UMAP are needed for alignment, so the AnnData is opened in backed mode and the count matrix stays on disk.

In [4]:

adata = ad.read_h5ad(H5AD_PATH, backed="r")
obs = adata.obs.copy()
spatial_um = np.asarray(adata.obsm["spatial"])
umap = np.asarray(adata.obsm["X_umap"])
colors = dict(zip(obs[CLUSTER_KEY].cat.categories, adata.uns[f"{CLUSTER_KEY}_colors"], strict=True))
adata.file.close()
print(f"{len(obs):,} cells, {obs['slice'].nunique()} slices, {adata.n_vars} genes")

adata = ad.read_h5ad(H5AD_PATH, backed="r") obs = adata.obs.copy() spatial_um = np.asarray(adata.obsm["spatial"]) umap = np.asarray(adata.obsm["X_umap"]) colors = dict(zip(obs[CLUSTER_KEY].cat.categories, adata.uns[f"{CLUSTER_KEY}_colors"], strict=True)) adata.file.close() print(f"{len(obs):,} cells, {obs['slice'].nunique()} slices, {adata.n_vars} genes")

3,828,568 cells, 44 slices, 480 genes

2. Leiden clusters¶

We use the leiden_1.0 clusters. Because the clusters are shared across slices, the same cluster appears in neighboring slices and can serve as a landmark for alignment.

In [5]:

print(f"{obs[CLUSTER_KEY].nunique()} Leiden clusters ({CLUSTER_KEY})")
rng = np.random.default_rng(0)
sub = rng.choice(len(obs), 300_000, replace=False)
fig, axes = plt.subplots(1, 2, figsize=(15, 6))
axes[0].scatter(
    umap[sub, 0],
    umap[sub, 1],
    c=obs[CLUSTER_KEY].iloc[sub].map(colors),
    s=0.3,
    linewidths=0,
    rasterized=True,
)
axes[0].set(title=f"UMAP, {CLUSTER_KEY}", xticks=[], yticks=[])
n_per_slice = obs.groupby("slice").size()
axes[1].bar(n_per_slice.index, n_per_slice.values, color="#4c78a8")
axes[1].set(xlabel="slice (anterior -> posterior)", ylabel="cells", title="cells per slice")
plt.tight_layout()
plt.show()

print(f"{obs[CLUSTER_KEY].nunique()} Leiden clusters ({CLUSTER_KEY})") rng = np.random.default_rng(0) sub = rng.choice(len(obs), 300_000, replace=False) fig, axes = plt.subplots(1, 2, figsize=(15, 6)) axes[0].scatter( umap[sub, 0], umap[sub, 1], c=obs[CLUSTER_KEY].iloc[sub].map(colors), s=0.3, linewidths=0, rasterized=True, ) axes[0].set(title=f"UMAP, {CLUSTER_KEY}", xticks=[], yticks=[]) n_per_slice = obs.groupby("slice").size() axes[1].bar(n_per_slice.index, n_per_slice.values, color="#4c78a8") axes[1].set(xlabel="slice (anterior -> posterior)", ylabel="cells", title="cells per slice") plt.tight_layout() plt.show()

38 Leiden clusters (leiden_1.0)

3. Select the cells to align¶

Two sets of cells are left out:

• Slice 44, a small fragment at the end of the series that barely overlaps its neighbor and has poor clusters. Its position in the stack is left empty, so the remaining slices keep their true spacing.
• A second, separate piece of tissue on slice 13 (flagged in extra_tissue_piece). Only the main piece is aligned.

In [6]:

keep = (~obs["slice"].isin(DROP_SLICES) & ~obs["extra_tissue_piece"]).to_numpy()

cells = ad.AnnData(obs=obs.loc[keep, [CLUSTER_KEY, "slice"]].copy())
cells.obsm["spatial"] = spatial_um[keep]  # Xenium positions are already in um
slices = cells.obs["slice"].to_numpy()
slice_ids = sorted(np.unique(slices))
print(f"{cells.n_obs:,} of {len(obs):,} cells kept, {len(slice_ids)} slices")

keep = (~obs["slice"].isin(DROP_SLICES) & ~obs["extra_tissue_piece"]).to_numpy() cells = ad.AnnData(obs=obs.loc[keep, [CLUSTER_KEY, "slice"]].copy()) cells.obsm["spatial"] = spatial_um[keep] # Xenium positions are already in um slices = cells.obs["slice"].to_numpy() slice_ids = sorted(np.unique(slices)) print(f"{cells.n_obs:,} of {len(obs):,} cells kept, {len(slice_ids)} slices")

3,787,366 of 3,828,568 cells kept, 43 slices

4. Align, step 1: cluster centroids¶

dega.align.calc_landmarks places one landmark per (slice, cluster) at the cluster's centroid. We also add a tissue-center landmark per slice. dega.align.calc_alignment_transform then fits a rigid Procrustes transform for each slice to its neighbors (within alignment_window slices), chained out from the middle slice. A cluster's centroid in one slice should land on top of the same cluster's centroid in the next.

In [7]:

landmarks = dega.align.calc_landmarks(cells, cluster_key=CLUSTER_KEY, slice_attr="slice")

rows = []
for s in slice_ids:
    xy = cells.obsm["spatial"][slices == s]
    (x_lo, x_hi), (y_lo, y_hi) = np.percentile(xy[:, 0], [1, 99]), np.percentile(xy[:, 1], [1, 99])
    rows.append(
        {
            "slice": s,
            "label": "tissue_center",
            "x": (x_lo + x_hi) / 2,
            "y": (y_lo + y_hi) / 2,
            "count": np.nan,
            "source": "manual",
        }
    )
landmarks = pd.concat([landmarks, pd.DataFrame(rows)], ignore_index=True)

procrustes = dega.align.calc_alignment_transform(
    landmarks,
    slice_attr="slice",
    reference=len(slice_ids) // 2,
    method="procrustes",
    **PROCRUSTES_KW,
)
print(f"{len(landmarks):,} landmarks, reference slice index {procrustes.reference}")

landmarks = dega.align.calc_landmarks(cells, cluster_key=CLUSTER_KEY, slice_attr="slice") rows = [] for s in slice_ids: xy = cells.obsm["spatial"][slices == s] (x_lo, x_hi), (y_lo, y_hi) = np.percentile(xy[:, 0], [1, 99]), np.percentile(xy[:, 1], [1, 99]) rows.append( { "slice": s, "label": "tissue_center", "x": (x_lo + x_hi) / 2, "y": (y_lo + y_hi) / 2, "count": np.nan, "source": "manual", } ) landmarks = pd.concat([landmarks, pd.DataFrame(rows)], ignore_index=True) procrustes = dega.align.calc_alignment_transform( landmarks, slice_attr="slice", reference=len(slice_ids) // 2, method="procrustes", **PROCRUSTES_KW, ) print(f"{len(landmarks):,} landmarks, reference slice index {procrustes.reference}")

1,278 landmarks, reference slice index 22

5. Align, step 2: neighborhoods¶

A neighborhood here is the alpha shape of all cells of one cluster in one slice, computed in each slice's native frame. dega.align.neighborhood_alignment starts from the Procrustes pose and nudges each slice (a bounded rotation and translation) to maximize the overlap of matching neighborhoods in nearby slices. Centroids average away a region's shape and orientation; overlapping the shapes recovers it.

On a laptop the alpha shapes take a few minutes and the refinement about 15-20 minutes.

In [8]:

%%time
shapes = dega.nbhd.alpha_shape_cell_clusters_by_slice(
    cells,
    cluster_attr=CLUSTER_KEY,
    slice_attr="slice",
    z_attr=None,
    alphas=(ALPHA_UM,),
    progress_every=10,
)
print(f"{len(shapes)} neighborhoods across {shapes['slice_id'].nunique()} slices")

%%time shapes = dega.nbhd.alpha_shape_cell_clusters_by_slice( cells, cluster_attr=CLUSTER_KEY, slice_attr="slice", z_attr=None, alphas=(ALPHA_UM,), progress_every=10, ) print(f"{len(shapes)} neighborhoods across {shapes['slice_id'].nunique()} slices")

alpha shapes: slice 10/43 (10)

alpha shapes: slice 20/43 (28)

alpha shapes: slice 30/43 (38)

alpha shapes: slice 40/43 (19)

alpha shapes: done, 976 neighborhoods across 43 slices
976 neighborhoods across 43 slices
CPU times: user 1min 21s, sys: 1.59 s, total: 1min 22s
Wall time: 1min 23s

In [9]:

%%time
neighborhood = dega.align.neighborhood_alignment(
    shapes,
    procrustes,
    slice_attr="slice_id",
    cluster_attr="cluster_id",
    progress_every=20,
    **NBHD_KW,
)

%%time neighborhood = dega.align.neighborhood_alignment( shapes, procrustes, slice_attr="slice_id", cluster_attr="cluster_id", progress_every=20, **NBHD_KW, )

[neighborhood_alignment +   0.0s] phase 1: chain-walk over 43 slices (window=4)

[neighborhood_alignment +  36.2s] 20 slice-optimizations

[neighborhood_alignment +  98.5s] 40 slice-optimizations

[neighborhood_alignment + 104.3s] phase 2: sweep 1/4

[neighborhood_alignment + 125.3s] 60 slice-optimizations

[neighborhood_alignment + 179.2s] 80 slice-optimizations

[neighborhood_alignment + 218.6s] 100 slice-optimizations

[neighborhood_alignment + 260.6s] 120 slice-optimizations

[neighborhood_alignment + 269.0s] phase 2: sweep 2/4

[neighborhood_alignment + 284.9s] 140 slice-optimizations

[neighborhood_alignment + 337.4s] 160 slice-optimizations

[neighborhood_alignment + 368.8s] 180 slice-optimizations

[neighborhood_alignment + 416.3s] 200 slice-optimizations

[neighborhood_alignment + 429.5s] phase 2: sweep 3/4

[neighborhood_alignment + 440.4s] 220 slice-optimizations

[neighborhood_alignment + 482.7s] 240 slice-optimizations

[neighborhood_alignment + 515.1s] 260 slice-optimizations

[neighborhood_alignment + 573.3s] 280 slice-optimizations

[neighborhood_alignment + 591.4s] phase 2: sweep 4/4

[neighborhood_alignment + 597.8s] 300 slice-optimizations

[neighborhood_alignment + 635.2s] 320 slice-optimizations

[neighborhood_alignment + 672.8s] 340 slice-optimizations

[neighborhood_alignment + 723.3s] 360 slice-optimizations

[neighborhood_alignment + 748.6s] optimization complete

CPU times: user 12min 11s, sys: 5.14 s, total: 12min 16s
Wall time: 12min 32s

6. Compare the two steps¶

For each pair of adjacent slices we compute the intersection-over-union (IoU) of each cluster's alpha shapes after alignment, averaged over the clusters the two slices share. Higher means matching regions overlap better.

In [10]:

# {slice: {cluster: 2D polygon}} in each slice's native frame
native = {}
for (s, cl), g in shapes.groupby(["slice_id", "cluster_id"]):
    native.setdefault(s, {})[str(cl)] = shapely.make_valid(
        shapely.force_2d(shapely.union_all(g.geometry.to_numpy()))
    )

alignments = {"centroids": procrustes.transforms, "neighborhood": neighborhood.transforms}


def adjacent_iou(transforms):
    out = []
    for a, b in itertools.pairwise(slice_ids):
        ious = []
        for cl in set(native[a]) & set(native[b]):
            ga = shapely.transform(native[a][cl], transforms[a].apply)
            gb = shapely.transform(native[b][cl], transforms[b].apply)
            inter = ga.intersection(gb).area
            union = ga.area + gb.area - inter
            ious.append(inter / union if union > 0 else 0.0)
        out.append(np.mean(ious) if ious else np.nan)
    return out


iou = pd.DataFrame({"slice": slice_ids[:-1]})
for name, t in alignments.items():
    iou[name] = adjacent_iou(t)

fig, ax = plt.subplots(figsize=(11, 3.5))
ax.plot(iou["slice"], iou["centroids"], "o-", ms=3, label="step 1: centroids")
ax.plot(iou["slice"], iou["neighborhood"], "o-", ms=3, label="step 2: + neighborhoods")
ax.set(xlabel="slice (IoU with the next slice)", ylabel="mean cluster IoU", ylim=(0, 1))
ax.legend(frameon=False)
plt.tight_layout()
plt.show()
iou[["centroids", "neighborhood"]].mean().round(3)

# {slice: {cluster: 2D polygon}} in each slice's native frame native = {} for (s, cl), g in shapes.groupby(["slice_id", "cluster_id"]): native.setdefault(s, {})[str(cl)] = shapely.make_valid( shapely.force_2d(shapely.union_all(g.geometry.to_numpy())) ) alignments = {"centroids": procrustes.transforms, "neighborhood": neighborhood.transforms} def adjacent_iou(transforms): out = [] for a, b in itertools.pairwise(slice_ids): ious = [] for cl in set(native[a]) & set(native[b]): ga = shapely.transform(native[a][cl], transforms[a].apply) gb = shapely.transform(native[b][cl], transforms[b].apply) inter = ga.intersection(gb).area union = ga.area + gb.area - inter ious.append(inter / union if union > 0 else 0.0) out.append(np.mean(ious) if ious else np.nan) return out iou = pd.DataFrame({"slice": slice_ids[:-1]}) for name, t in alignments.items(): iou[name] = adjacent_iou(t) fig, ax = plt.subplots(figsize=(11, 3.5)) ax.plot(iou["slice"], iou["centroids"], "o-", ms=3, label="step 1: centroids") ax.plot(iou["slice"], iou["neighborhood"], "o-", ms=3, label="step 2: + neighborhoods") ax.set(xlabel="slice (IoU with the next slice)", ylabel="mean cluster IoU", ylim=(0, 1)) ax.legend(frameon=False) plt.tight_layout() plt.show() iou[["centroids", "neighborhood"]].mean().round(3)

Out[10]:

centroids       0.336
neighborhood    0.426
dtype: float64

7. Apply both alignments to every cell¶

Each slice's transform is applied to its cells. The position along the stack is (slice - reference) * 100 µm, so the dropped slice leaves a gap instead of shifting the later slices.

The plots show all slices stacked on top of each other (left) and from the side (right). After the neighborhood step the clusters line up more tightly from slice to slice.

In [11]:

xy_native = cells.obsm["spatial"]
stack_um = (slices - procrustes.reference) * Z_SPACING_UM

aligned_xy = {}
for name, transforms in alignments.items():
    xy = np.empty_like(xy_native)
    for s in slice_ids:
        m = slices == s
        xy[m] = transforms[s].apply(xy_native[m])
    aligned_xy[name] = xy

point_colors = cells.obs[CLUSTER_KEY].map(colors).to_numpy()
fig, axes = plt.subplots(2, 2, figsize=(16, 12))
for row, (name, xy) in enumerate(aligned_xy.items()):
    axes[row, 0].scatter(xy[:, 0], xy[:, 1], c=point_colors, s=0.02, linewidths=0, rasterized=True)
    axes[row, 0].set(title=f"{name}: all slices, looking along the stack", aspect="equal")
    axes[row, 0].invert_yaxis()
    axes[row, 1].scatter(stack_um, xy[:, 1], c=point_colors, s=0.02, linewidths=0, rasterized=True)
    axes[row, 1].set(title=f"{name}: side view", xlabel="stack (um)")
    axes[row, 1].invert_yaxis()
plt.tight_layout()
plt.show()

xy_native = cells.obsm["spatial"] stack_um = (slices - procrustes.reference) * Z_SPACING_UM aligned_xy = {} for name, transforms in alignments.items(): xy = np.empty_like(xy_native) for s in slice_ids: m = slices == s xy[m] = transforms[s].apply(xy_native[m]) aligned_xy[name] = xy point_colors = cells.obs[CLUSTER_KEY].map(colors).to_numpy() fig, axes = plt.subplots(2, 2, figsize=(16, 12)) for row, (name, xy) in enumerate(aligned_xy.items()): axes[row, 0].scatter(xy[:, 0], xy[:, 1], c=point_colors, s=0.02, linewidths=0, rasterized=True) axes[row, 0].set(title=f"{name}: all slices, looking along the stack", aspect="equal") axes[row, 0].invert_yaxis() axes[row, 1].scatter(stack_um, xy[:, 1], c=point_colors, s=0.02, linewidths=0, rasterized=True) axes[row, 1].set(title=f"{name}: side view", xlabel="stack (um)") axes[row, 1].invert_yaxis() plt.tight_layout() plt.show()

8. Write the CellCloud DegaFiles¶

The DegaFiles for this dataset are already hosted on GitHub at broadinstitute/celldega_xenium_E15.5_cell_cloud, so you don't need to write them to follow the Viz notebook. The code below shows how they were made.

CellCloud draws every cell as a point in 3D. Its orbit camera spins around the vertical (data-Z) axis, so side_view_coords puts the stack along data-X (anterior on the left) and the dorsal-ventral axis along Z (dorsal up). Points are drawn in file order without depth sorting, so the cells are written anterior first and posterior last, to suit the Viz notebook's camera.

dega.align.write_alignment_point_cloud writes both alignments as named variants of one directory. The first call creates the directory, including the cluster colors and the counts of all 480 genes (cbg/, one file per gene). The second call only adds the centroids positions.

def side_view_coords(xy_um, stack_um):
    """(x, y, stack) -> (stack, x, y), centered on the origin, with dorsal up."""
    pts = np.column_stack([xy_um, stack_um])
    center = (pts.min(axis=0) + pts.max(axis=0)) / 2
    x, y, s = (pts - center).T
    return np.column_stack([s, x, y])


CELL_CLOUD_DIR = Path("DegaFiles/xenium_E15.5_cell_cloud")
draw_order = np.argsort(slices, kind="stable")  # posterior slices last

cloud = ad.read_h5ad(H5AD_PATH)[keep][draw_order].copy()  # raw counts in X
cloud.obs = cloud.obs[[CLUSTER_KEY]]
for name in ["neighborhood", "centroids"]:
    # rounded to 1 um: sub-micron precision isn't needed and the files shrink by ~30%
    xyz = np.round(side_view_coords(aligned_xy[name][draw_order], stack_um[draw_order]))
    cloud.obsm["spatial"], cloud.obs["Z"] = xyz[:, :2], xyz[:, 2]
    dega.align.write_alignment_point_cloud(
        cloud, CELL_CLOUD_DIR, alignment_name=name, z_key="Z", cluster_key=CLUSTER_KEY, overwrite=True
    )

Next¶

Broad-Allen_Dev-Mouse-Brain_E15.5_Viz loads these DegaFiles into a 3D CellCloud linked to the cluster composition of every slice and to a Clustergram of gene expression.