Quantitative Assessment of Mesoscale Cellular Order and Organization in the Mouse Hippocampus.

Code ↔ Paper

3 matches between paragraphs of the paper and lines of its authors' code, computed by the harvester (lexical-v1). Click a colored paragraph or line to see its counterpart.

The 3 matches

[1] § Results › Intensity Processing and Segmentation of Nuclei in High‐Density Brain Areas ↔ code/FunctionsClean.py, lines 68 to 107 · score 0.60 · background illumination, Gaussian filter, dilation, fitting · Report a mismatch
[2] § Material and Methods › Image Acquisition ↔ code/cbhs3d_onefile.py, lines 95 to 96 · score 0.54 · Zeiss LSM, 0.1135 um · Report a mismatch
[3] § Results › Intensity Processing and Segmentation of Nuclei in High‐Density Brain Areas ↔ code/cbhs3d_onefile.py, lines 357 to 395 · score 0.50 · Gaussian filter, noise, signals, dilation, erosion, components · Report a mismatch

Paper

Loaded from Europe PMC by your browser; OSCR does not store it: doi.org · Europe PMC

The paper is loaded when this pane is shown.

The paper is loaded by your browser, with JavaScript: read it at doi.org or on Europe PMC.

The authors' code

Python · 1,250 lines · 55 KB · no license · 2 matches

#!/usr/bin/env python
r"""CBHS 3D pipeline, single file: raw TIFF z-stack -> 3D segmented nuclear masks + shape metrics.
Implements the README "Current pipeline" (the chain that produced ``M_synth3_yuki_*`` and
``M_synth3_yuki_anchored``) and the README "Shape metrics" section, with no imports from
the ``cbhs3d`` package -- every function is inlined here.
raw z-slices --(1) cellpose 2D, independent-----> M_indep.tiff
M_indep + raw --(2) synthetic DAPI v3-----------> synth_dapi_v3.tiff
synth --(3) StarDist 3D (own env)-------> M_stardist_v3.tiff
M_stardist --(4) Yukimatsu refine (reduced)--> M_synth3_yuki_reduced.tiff
+ M_indep --(5) anchor 3D->2D + finish------> M_synth3_yuki_reduced_anchored.tiff
final masks --(6) shape metrics---------------> per_nucleus.csv, per_slice.csv
Stage 3 needs TensorFlow, which the analysis environment does not have. Pass
``--stardist-python`` and this file re-invokes *itself* with that interpreter for stage 3
only (``--stardist-stage``); omit it and stage 3 runs in-process instead. Stages are
cached: an existing output file skips its stage unless ``--force``.
PROVENANCE (original script -> section below)
Sec. 2 Cellpose front-ends ............ Cellpose3DSegment.py, Compare3DvsMean2D.cellpose_segment
Sec. 5a round_caps, split_touching_shape_aware, relabel_compact_1n
Yukimatsu.py
Sec. 5b cleanup_labels ................. YukimatsuOnSynth2Stardist.cleanup_morphology
VARIANTS ....................... YukimatsuOnSynth3.py
Sec. 5c link_2d_to_3d .................. HybridSegment.build_hybrid
--- included with explicit approval, not from the four scripts above ---
Sec. 3 synthetic DAPI ................. BuildSyntheticDapi.py / _v2.py / _v3.py, PipelineV2.py st.1
Sec. 4 run_stardist3d ................. Stardist3DSegment.py, Stardist3DOnSynth*.py, PipelineV2 st.2
Sec. 5d anchor_to_2d, anchor_and_finish AnchorTo2D.py
Sec. 6 metrics_3d, metrics_2d_* ....... Compare3DvsMean2D.py
per_slice_table ................ SliceSampling3Dvs2D.build_per_slice
Sec. 1 z-stack loading, MIPs, save .... library glue (no original-script equivalent)
NOT included: audit_fragmentation (AuditHybridYuki.py) -- not needed by this chain.
Usage
python cbhs3d_onefile.py --raw-dir "<...IHCFinal_6wETHO_4_1 c=0.tiff_files>" \
--out-dir out \
--cellpose-weights <cpsam_20250704_112416> \
--stardist-model-dir <Ferme-et-al_2025/stardist> \
--stardist-python C:\Users\Admin\.conda\envs\stardist-gpu\python.exe \
--synth v3 --variant reduced --anchor --metrics
Environments
analysis: numpy, scipy, scikit-image, pandas, tifffile, tqdm, cellpose, torch (GPU rounding)
stardist-gpu: Python 3.10, TF 2.10, CUDA 11.2, stardist, csbdeep (stage 3 only)
"""
import argparse
import glob
import json
import math
import os
import re
import subprocess
import sys
os.environ.setdefault("KMP_DUPLICATE_LIB_OK", "TRUE")
import numpy as np
import tifffile
from scipy.ndimage import (binary_dilation, binary_erosion, binary_fill_holes,
distance_transform_edt, find_objects, gaussian_filter,
label as cc_label, zoom)
from scipy.spatial import ConvexHull
from skimage.feature import peak_local_max
from skimage.measure import marching_cubes, mesh_surface_area, regionprops
from skimage.segmentation import watershed
from tqdm import tqdm
try: # only Sec. 6 needs pandas; the stardist env may lack it
import pandas as pd
except ImportError:
pd = None
try: # re-exported at scipy.spatial only from SciPy 1.8
from scipy.spatial import QhullError
except ImportError:
from scipy.spatial.qhull import QhullError
try: # torch missing or broken -> CPU morphology path
import torch
import torch.nn.functional as F
_HAS_CUDA = torch.cuda.is_available()
except Exception:
torch = None
F = None
_HAS_CUDA = False
# ============================================================================
# 0. CONFIG -- calibration and parameter presets (config.py)
# Every preset is a plain dict; copy and edit rather than mutating in place.
# ============================================================================
# (Z, Y, X) voxel size in micrometres -- Zeiss LSM, 40x, 0.33 um z-step
VOXEL_UM = (0.33, 0.1135, 0.1135)
# Cellpose 2D (custom cpsam model)
CELLPOSE = dict(diameter=45, flow_threshold=0.6, min_size=30, stitch_threshold=0.4)
# StarDist 3D (Ferme et al. 2025, Zenodo 10.5281/zenodo.15120155, model_A)
STARDIST = dict(model_name="model_A", n_tiles=(2, 8, 8),
prob_thresh=None, nms_thresh=None, normalize_pct=(1.0, 99.8))
# Synthetic DAPI (2D masks -> 3D-segmentable volume)
SYNTH_V1 = dict(kind="binary", intensity_interior=255, intensity_rim=160, rim_width_vox=2)
SYNTH_V2 = dict(kind="masked", xy_erosion_vox=2, z_dilation=True, xy_smooth_sigma=0.7,
pct_lo=2.0, pct_hi=98.0, noise_sigma=5.0, seed=42)
SYNTH_V3 = dict(SYNTH_V2, xy_erosion_vox=1) # reduced-rim variant (current)
# Yukimatsu cap rounding
ROUNDING_DEFAULT = dict(
z_boost=2.6, max_z_scale=4, gaussian_sigma=1.1,
close_radius_xy_um=0.35, close_radius_z_um=0.90, open_radius_um=0.30,
virtual_z_margin_um=1.2, imagej_interp=True, z_push_um=0.30,
)
# reduced rounding: no z-push, smaller Gaussian and z-closing radii, no opening (synth v3 variant)
ROUNDING_REDUCED = dict(ROUNDING_DEFAULT, z_push_um=0.0, gaussian_sigma=0.6,
close_radius_xy_um=0.35, close_radius_z_um=0.45,
open_radius_um=0.0, virtual_z_margin_um=0.4)
# rounding applied after anchoring to the 2D masks; closes residual slice steps
ROUNDING_POST_ANCHOR = dict(ROUNDING_DEFAULT, z_push_um=0.25, gaussian_sigma=0.6,
close_radius_xy_um=0.20, close_radius_z_um=0.70,
open_radius_um=0.0, virtual_z_margin_um=1.0)
# rounding with larger radii, applied to Cellpose-stitched input (Compare3DvsMean2D M_final)
ROUNDING_STITCHED = dict(ROUNDING_DEFAULT, gaussian_sigma=2.0,
close_radius_xy_um=0.80, close_radius_z_um=1.50,
open_radius_um=0.70, virtual_z_margin_um=1.5, z_push_um=0.45)
# Yukimatsu shape-aware split (intensity-free, DT markers + neck/shape guards)
SPLIT_DEFAULT = dict(min_peak_sep_um=1.0, max_peaks=5, min_rel_dt_at_peak=0.22,
neck_rel=0.35, min_child_voxels=2600, min_child_vol_ratio=0.25,
elong_merge_if_ratio_below=1.15)
SPLIT_STITCHED = dict(min_peak_sep_um=0.8, max_peaks=8, min_rel_dt_at_peak=0.18,
neck_rel=0.25, min_child_voxels=1500, min_child_vol_ratio=0.18,
elong_merge_if_ratio_below=1.10)
# Per-label cleanup (largest CC + hole fill + ellipsoid opening + size filter)
CLEANUP_DEFAULT = dict(open_radius_um=1.2, min_voxels=200)
CLEANUP_REDUCED = dict(open_radius_um=0.5, min_voxels=200)
CLEANUP_ANCHOR = dict(open_radius_um=0.6, min_voxels=200) # removes shelves 1-3 slices thick
CLEANUP_NO_OPENING = dict(open_radius_um=0.0, min_voxels=200)
# Named end-to-end variants (name -> rounding, cleanup); None = no rounding.
# Former names: default = published, reduced = gentle, no_rounding = minimal.
VARIANTS = {
"default": (ROUNDING_DEFAULT, CLEANUP_DEFAULT),
"reduced": (ROUNDING_REDUCED, CLEANUP_REDUCED),
"no_rounding": (None, CLEANUP_NO_OPENING),
}
# 2D <-> 3D linking
HYBRID_OVERLAP_THR = 0.30 # Cellpose blob must have >=30 % of its pixels in one 3D label
ANCHOR_PAD_PX = 20 # Cellpose blob may exceed the 3D bbox by this many pixels
# ============================================================================
# 1. I/O (io.py -- library glue)
# ============================================================================
Z_INDEX_PATTERN = re.compile(r"b0z(\d+)c", re.IGNORECASE)
def load_zstack_from_folder(folder, pattern=Z_INDEX_PATTERN):
"""Stack per-z TIFFs (Zeiss export ``..._b0z<N>c0...tiff``) into a (Z, Y, X) volume.
RGB/RGBA exports with R==G==B are reduced to channel 0."""
files = [f for f in glob.glob(os.path.join(folder, "*.tif*"))
if "_metadata" not in os.path.basename(f)]
indexed = []
for f in files:
m = pattern.search(os.path.basename(f))
if m:
indexed.append((int(m.group(1)), f))
if not indexed:
raise FileNotFoundError(f"no per-z tiffs matching {pattern.pattern!r} in {folder}")
indexed.sort()
slices = []
for _, f in indexed:
a = tifffile.imread(f)
if a.ndim == 3 and a.shape[-1] in (3, 4):
a = a[..., 0]
slices.append(a)
return np.stack(slices, axis=0)
def imread_3d(path):
arr = tifffile.imread(path)
if arr.ndim != 3:
raise ValueError(f"{os.path.basename(path)} is not 3D (ndim={arr.ndim})")
return arr
def load_raw(path):
"""A folder of per-z TIFFs or a single 3D TIFF."""
return load_zstack_from_folder(path) if os.path.isdir(path) else imread_3d(path)
def n_labels(masks):
return int(np.unique(masks).size - (1 if (masks == 0).any() else 0))
def relabel_compact_1n(masks, dtype=np.uint32):
"""Relabel to 1..N preserving label order; keeps all labels. [Yukimatsu.relabel_compact_1N]"""
labels = np.unique(masks)
labels = labels[labels > 0]
if labels.size == 0:
return masks.astype(dtype, copy=False)
lut = np.zeros(int(masks.max()) + 1, dtype=dtype)
lut[labels] = np.arange(1, labels.size + 1, dtype=dtype)
return lut[masks]
def ensure_unique(path):
if not os.path.exists(path):
return path
stem, ext = os.path.splitext(path)
k = 2
while os.path.exists(f"{stem}__v{k}{ext}"):
k += 1
return f"{stem}__v{k}{ext}"
def save_mips(stem, vol):
"""Write XY / XZ / YZ maximum-intensity projections next to ``stem``."""
for name, axis in (("XY", 0), ("XZ", 1), ("YZ", 2)):
tifffile.imwrite(ensure_unique(f"{stem}_MIP_{name}.tiff"),
np.max(vol, axis=axis).astype(np.uint16, copy=False))
def save_labels(path, masks):
tifffile.imwrite(path, np.asarray(masks).astype(np.uint32, copy=False))
# ============================================================================
# 2. CELLPOSE FRONT-ENDS [Cellpose3DSegment.py, Compare3DvsMean2D.cellpose_segment]
# Require the ``cellpose`` package (CUDA optional).
# ============================================================================
def _cellpose_model(weights_path, gpu=True):
from cellpose import models
return models.CellposeModel(gpu=gpu, pretrained_model=weights_path)
def cellpose_segment_independent(stack, weights_path, diameter=CELLPOSE["diameter"],
flow_threshold=CELLPOSE["flow_threshold"],
min_size=CELLPOSE["min_size"], gpu=True):
"""2D per slice, no stitching. Labels are offset per slice so they are globally unique.
This is stage 1 of the current pipeline."""
model = _cellpose_model(weights_path, gpu)
out = model.eval([stack[z] for z in range(stack.shape[0])], diameter=diameter,
do_3D=False, flow_threshold=flow_threshold, min_size=min_size)
per_slice = out[0] if isinstance(out, tuple) else out
full = np.zeros(stack.shape, dtype=np.uint32)
offset = 0
for z, m in enumerate(per_slice):
m = np.asarray(m, dtype=np.uint32)
nz = m > 0
if nz.any():
m[nz] += offset
offset = int(m.max())
full[z] = m
return full
def cellpose_segment_stitched(stack, weights_path, stitch_threshold=CELLPOSE["stitch_threshold"],
diameter=CELLPOSE["diameter"], flow_threshold=CELLPOSE["flow_threshold"],
min_size=CELLPOSE["min_size"], gpu=True):
"""2D per slice + z-stitching -> labels consistent across z. Returns uint32 (Z, Y, X).
This is the legacy (paper Fig. 3 F/G) stage 1."""
model = _cellpose_model(weights_path, gpu)
out = model.eval(stack, z_axis=0, diameter=diameter, do_3D=False,
stitch_threshold=stitch_threshold, flow_threshold=flow_threshold,
min_size=min_size)
masks = out[0] if isinstance(out, tuple) else out
return np.asarray(masks).astype(np.uint32)
def cellpose_segment_3d(stack, weights_path, voxel_um=VOXEL_UM, diameter=CELLPOSE["diameter"],
flow_threshold=CELLPOSE["flow_threshold"], min_size=CELLPOSE["min_size"],
gpu=True):
"""Native Cellpose 3D (do_3D=True): XY/XZ/YZ flows combined. Slow (30-60 min on GPU)."""
model = _cellpose_model(weights_path, gpu)
anisotropy = voxel_um[0] / voxel_um[1]
out = model.eval(stack, z_axis=0, diameter=diameter, do_3D=True,
flow_threshold=flow_threshold, min_size=min_size, anisotropy=anisotropy)
masks = out[0] if isinstance(out, tuple) else out
return np.asarray(masks).astype(np.uint32)
CELLPOSE_MODES = {
"independent": cellpose_segment_independent, # current pipeline, stage 1
"stitched": cellpose_segment_stitched, # legacy chain
"3d": cellpose_segment_3d,
}
# ============================================================================
# 3. SYNTHETIC DAPI [BuildSyntheticDapi.py / _v2.py / _v3.py, PipelineV2.py stage 1]
# Volumes built from the independent 2D masks, so that a 3D segmenter (StarDist 3D)
# can link slices while the in-plane shape stays that of the 2D masks.
#
# v1 binary interior + darker rim, hard edges -> build_synth_dapi_binary
# v2 raw DAPI masked by 2-vox-eroded blobs + refinements -> build_synth_dapi_masked(2)
# v3 as v2 but 1-vox erosion (reduced rim, current) -> build_synth_dapi_masked(1)
# ============================================================================
def z_dilate_labels(M):
"""Propagate every 2D blob one slice up and down, into empty voxels only.
Recovers faint cap slices without bridging stacked nuclei."""
Z = M.shape[0]
out = M.copy()
for z in range(Z - 1):
fill = (out[z + 1] == 0) & (M[z] > 0)
out[z + 1] = np.where(fill, M[z], out[z + 1])
for z in range(Z - 1, 0, -1):
fill = (out[z - 1] == 0) & (M[z] > 0)
out[z - 1] = np.where(fill, M[z], out[z - 1])
return out
def erode_blobs_xy(M, iterations):
"""Per-blob 2D erosion; returns a boolean foreground with gaps between neighbours."""
Z = M.shape[0]
fg = np.zeros(M.shape, dtype=bool)
for z in tqdm(range(Z), desc="erode xy", unit="slice", leave=False):
sl2 = M[z]
for idx, sl in enumerate(find_objects(sl2)):
if sl is None:
continue
m = sl2[sl] == idx + 1
if m.any():
fg[z][sl] |= binary_erosion(m, iterations=iterations) if iterations > 0 else m
return fg
def build_synth_dapi_binary(M_indep, intensity_interior=SYNTH_V1["intensity_interior"],
intensity_rim=SYNTH_V1["intensity_rim"],
rim_width_vox=SYNTH_V1["rim_width_vox"]):
"""v1: interior=255, rim=160; the rim erosion doubles as the gap to neighbours."""
Z = M_indep.shape[0]
synth = np.zeros(M_indep.shape, dtype=np.uint8)
for z in tqdm(range(Z), desc="synth v1", unit="slice", leave=False):
sl2 = M_indep[z]
for idx, sl in enumerate(find_objects(sl2)):
if sl is None:
continue
m = sl2[sl] == idx + 1
if not m.any():
continue
interior = binary_erosion(m, iterations=rim_width_vox)
sub = synth[z][sl]
sub[m & ~interior] = intensity_rim
sub[interior] = intensity_interior
return synth
def build_synth_dapi_masked(raw, M_indep, xy_erosion_vox=SYNTH_V3["xy_erosion_vox"],
z_dilation=SYNTH_V2["z_dilation"],
xy_smooth_sigma=SYNTH_V2["xy_smooth_sigma"],
pct_lo=SYNTH_V2["pct_lo"], pct_hi=SYNTH_V2["pct_hi"],
noise_sigma=SYNTH_V2["noise_sigma"], seed=SYNTH_V2["seed"]):
"""v2/v3: raw DAPI inside eroded 2D blobs, per-3D-connected-component percentile
normalisation (preserves the natural z-cap fade), light XY smoothing and additive
Gaussian noise on signal voxels. Returns uint8."""
if raw.shape != M_indep.shape:
raise ValueError(f"shape mismatch raw {raw.shape} vs masks {M_indep.shape}")
M = z_dilate_labels(M_indep) if z_dilation else M_indep
fg = erode_blobs_xy(M, xy_erosion_vox)
synth = np.zeros(raw.shape, dtype=np.float32)
synth[fg] = raw[fg]
cc, _ = cc_label(fg, structure=np.ones((3, 3, 3)))
out = np.zeros_like(synth)
for idx, sl in enumerate(tqdm(find_objects(cc), desc="normalise cc", unit="cc", leave=False)):
if sl is None:
continue
m = cc[sl] == idx + 1
if not m.any():
continue
vals = synth[sl][m]
if vals.size < 5:
out[sl][m] = vals
continue
p_lo, p_hi = np.percentile(vals, [pct_lo, pct_hi])
if p_hi <= p_lo:
out[sl][m] = vals
continue
out[sl][m] = np.clip((vals - p_lo) / (p_hi - p_lo) * 240.0 + 15.0, 0, 255)
out = gaussian_filter(out, sigma=(0.0, xy_smooth_sigma, xy_smooth_sigma))
if noise_sigma > 0:
rng = np.random.default_rng(seed)
out = out + rng.normal(0.0, noise_sigma, out.shape).astype(np.float32) * (out > 0)
return np.clip(out, 0, 255).astype(np.uint8)
def build_synth_dapi(raw, M_indep, version="v3"):
"""Dispatch by version name: 'v1' | 'v2' | 'v3'."""
if version == "v1":
return build_synth_dapi_binary(M_indep, **{k: v for k, v in SYNTH_V1.items() if k != "kind"})
cfg = SYNTH_V2 if version == "v2" else SYNTH_V3
return build_synth_dapi_masked(raw, M_indep, **{k: v for k, v in cfg.items() if k != "kind"})
# ============================================================================
# 4. STARDIST 3D [Stardist3DSegment.py, Stardist3DOnSynth*.py, PipelineV2.py stage 2]
# Needs the ``stardist-gpu`` environment (Python 3.10, TF 2.10, CUDA 11.2,
# stardist, csbdeep). Nothing else in this file imports TensorFlow, and the
# imports below are function-local so the analysis env can load this module.
# ============================================================================
def run_stardist3d(stack, model_basedir, model_name=STARDIST["model_name"],
n_tiles=STARDIST["n_tiles"], prob_thresh=STARDIST["prob_thresh"],
nms_thresh=STARDIST["nms_thresh"], normalize_pct=STARDIST["normalize_pct"],
verbose=True):
"""Percentile-normalise, pad Z to a multiple of 4 (UNet alignment), predict, crop back.
Returns uint32 labels of the input shape."""
os.environ.setdefault("TF_CPP_MIN_LOG_LEVEL", "2")
os.environ.setdefault("TF_ENABLE_ONEDNN_OPTS", "0")
from csbdeep.utils import normalize
from stardist.models import StarDist3D
x = normalize(stack.astype(np.float32, copy=False), *normalize_pct, axis=None)
Z = x.shape[0]
Z_pad = ((Z + 3) // 4) * 4
if Z_pad != Z:
x = np.pad(x, ((0, Z_pad - Z), (0, 0), (0, 0)), mode="edge")
model = StarDist3D(None, name=model_name, basedir=model_basedir)
kw = dict(n_tiles=n_tiles, verbose=verbose)
if prob_thresh is not None:
kw["prob_thresh"] = prob_thresh
if nms_thresh is not None:
kw["nms_thresh"] = nms_thresh
labels, _ = model.predict_instances(x, **kw)
return np.asarray(labels[:Z]).astype(np.uint32)
# ============================================================================
# 5. YUKIMATSU -- label-volume refinement for anisotropic nuclear masks
#
# (a) round_caps z-upsample, z-push, Gaussian smooth, ellipsoidal
# close/open, max-pool back [Yukimatsu.py]
# split_touching_shape_aware DT watershed with neck / size / elongation guards
# [Yukimatsu.py]
# (b) cleanup_labels largest CC + hole fill + ellipsoid opening + size
# filter [YukimatsuOnSynth2Stardist.py]
# (c) link_2d_to_3d each 2D blob takes the label of the 3D object it
# overlaps most [HybridSegment.py]
# (d) anchor_to_2d back-match: rebuild every 3D label slice by slice
# from the 2D masks [AnchorTo2D.py]
# (e) refine chain: compact -> round -> [split] -> cleanup ->
# [anchor]
#
# All radii are in micrometres and converted with the (Z, Y, X) voxel size.
# ============================================================================
def _r(v, minimum=0):
return int(max(minimum, round(float(v))))
def ellipsoid_se(rz, ry, rx):
rz_s, ry_s, rx_s = max(rz, 1), max(ry, 1), max(rx, 1)
Z, Y, X = np.ogrid[-rz:rz + 1, -ry:ry + 1, -rx:rx + 1]
return (Z * Z) / (rz_s * rz_s) + (Y * Y) / (ry_s * ry_s) + (X * X) / (rx_s * rx_s) <= 1.0
def largest_cc(mask, connectivity=1):
structure = np.ones((3, 3, 3)) if connectivity == 3 else None
cc, n = cc_label(mask, structure=structure)
if n <= 1:
return mask
counts = np.bincount(cc.ravel())
counts[0] = 0
return cc == int(np.argmax(counts))
def _radii_vox(um, voxel_um, minimum=0):
return tuple(_r(um / v, minimum) for v in voxel_um)
# ---- GPU morphology (separable box approximation of the ellipsoid, as in the original) ----
def _gauss1d_gpu(s, dev):
if s <= 0:
return torch.tensor([1.0], device=dev)
r = max(1, int(3 * s + 0.5))
x = torch.arange(-r, r + 1, device=dev, dtype=torch.float32)
k = torch.exp(-(x * x) / (2 * s * s))
return k / k.sum()
def _blur3d_gpu(vol, s, dev):
k = _gauss1d_gpu(s, dev)
p = k.numel() // 2
vol = F.conv3d(F.pad(vol, (0, 0, 0, 0, p, p), mode="replicate"), k.view(1, 1, -1, 1, 1))
vol = F.conv3d(F.pad(vol, (0, 0, p, p, 0, 0), mode="replicate"), k.view(1, 1, 1, -1, 1))
vol = F.conv3d(F.pad(vol, (p, p, 0, 0, 0, 0), mode="replicate"), k.view(1, 1, 1, 1, -1))
return vol
def _morph_gpu(vol, rz, ry, rx, op, dev):
def pass1(v, r, axis):
if r <= 0:
return v
shape = [1, 1, 1, 1, 1]
shape[2 + axis] = 2 * r + 1
k = torch.ones(shape, device=dev)
pad = [0, 0, 0, 0, 0, 0]
pad[2 * (2 - axis)] = pad[2 * (2 - axis) + 1] = r
c = F.conv3d(F.pad(v, pad, mode="constant", value=0.0), k)
return (c > 0).float() if op == "dilate" else (c >= (2 * r + 1) - 1e-6).float()
v = pass1(vol, rz, 0)
v = pass1(v, ry, 1)
return pass1(v, rx, 2)
def _round_one_gpu(sub, z_scale, tri, rz_push, sigma, close_r, open_r, dev):
t = torch.from_numpy(sub.astype(np.float32))[None, None].to(dev, non_blocking=True)
if z_scale > 1:
t = F.interpolate(t, scale_factor=(z_scale, 1, 1),
mode="trilinear" if tri else "nearest",
align_corners=False if tri else None)
if rz_push > 0:
k = torch.ones((1, 1, 2 * rz_push + 1, 1, 1), device=dev)
t = (F.conv3d(F.pad(t, (0, 0, 0, 0, rz_push, rz_push)), k) > 0).float()
t = (_blur3d_gpu(t, sigma, dev) >= 0.5).float()
t = _morph_gpu(t, *close_r, "dilate", dev)
t = _morph_gpu(t, *close_r, "erode", dev)
if max(open_r) > 0:
t = _morph_gpu(t, *open_r, "erode", dev)
t = _morph_gpu(t, *open_r, "dilate", dev)
if z_scale > 1:
t = F.max_pool3d(t, kernel_size=(z_scale, 1, 1), stride=(z_scale, 1, 1))
return t.squeeze().detach().cpu().numpy() > 0.5
def _round_one_cpu(sub, z_scale, tri, rz_push, sigma, close_r, open_r):
if z_scale > 1:
hi = zoom(sub.astype(np.float32), zoom=(z_scale, 1, 1), order=1 if tri else 0) >= 0.5
else:
hi = sub.astype(bool, copy=False)
if rz_push > 0:
hi = binary_dilation(hi, structure=np.ones((2 * rz_push + 1, 1, 1), dtype=bool))
hi = gaussian_filter(hi.astype(np.float32), sigma=sigma, mode="nearest") >= 0.5
se = ellipsoid_se(*close_r)
hi = binary_erosion(binary_dilation(hi, structure=se), structure=se)
if max(open_r) > 0:
se2 = ellipsoid_se(*open_r)
hi = binary_dilation(binary_erosion(hi, structure=se2), structure=se2)
if z_scale <= 1:
return hi
Zhi = hi.shape[0]
out = np.zeros(sub.shape, dtype=bool)
for zi in range(sub.shape[0]):
z0 = max(0, min(Zhi, int(round(zi * z_scale))))
z1 = max(z0 + 1, min(Zhi, int(round((zi + 1) * z_scale))))
out[zi] = np.any(hi[z0:z1], axis=0)
return out
# ---- 5a. cap rounding [Yukimatsu.round_caps] ----
def round_caps(masks, voxel_um=VOXEL_UM, cfg=ROUNDING_DEFAULT, use_gpu=True,
xy_pad=4, progress=True):
"""Per-label z-cap rounding. Each label is upsampled along Z (``z_boost`` x anisotropy,
capped at ``max_z_scale``), optionally pushed along Z, Gaussian-smoothed, closed and
opened with ellipsoidal elements, then max-pooled back to the original z grid so the
rounded cap may spill into ``virtual_z_margin_um`` of neighbouring original slices.
Larger labels are written first; overlaps are first-come."""
masks = np.asarray(masks)
if masks.max() == 0:
return masks.copy()
z_um, y_um, x_um = map(float, voxel_um)
anis = max(1.0, z_um / min(y_um, x_um))
z_scale = int(min(cfg.get("max_z_scale", 4), max(1, round(anis * cfg.get("z_boost", 2.0)))))
z_hi = z_um / z_scale
sigma = float(cfg.get("gaussian_sigma", 1.8))
tri = bool(cfg.get("imagej_interp", True))
close_r = (_r(cfg.get("close_radius_z_um", 0.90) / z_hi, 1),
_r(cfg.get("close_radius_xy_um", 0.35) / y_um, 1),
_r(cfg.get("close_radius_xy_um", 0.35) / x_um, 1))
open_um = float(cfg.get("open_radius_um", 0.30))
open_r = (_r(open_um / z_hi), _r(open_um / y_um), _r(open_um / x_um))
rz_push = _r(cfg.get("z_push_um", 0.60) / z_hi)
spill_um = max(float(cfg.get("virtual_z_margin_um", 1.2)), z_um + 1e-6)
z_pad = max(1, int(np.ceil(spill_um / z_um)))
gpu = bool(use_gpu and _HAS_CUDA)
dev = torch.device("cuda") if gpu else None
out = np.zeros_like(masks)
sizes = np.bincount(masks.ravel())
bboxes = find_objects(masks)
labs = [i for i in range(1, sizes.size) if sizes[i] > 0 and bboxes[i - 1] is not None]
labs.sort(key=lambda L: sizes[L], reverse=True)
Z, Y, X = masks.shape
it = tqdm(labs, desc=f"round caps ({'GPU' if gpu else 'CPU'})", unit="obj", leave=False) \
if progress else labs
for lab in it:
sl = bboxes[lab - 1]
z0, z1 = max(0, sl[0].start - 1 - z_pad), min(Z, sl[0].stop + 1 + z_pad)
y0, y1 = max(0, sl[1].start - xy_pad), min(Y, sl[1].stop + xy_pad)
x0, x1 = max(0, sl[2].start - xy_pad), min(X, sl[2].stop + xy_pad)
slc = (slice(z0, z1), slice(y0, y1), slice(x0, x1))
sub = masks[slc] == lab
if not sub.any():
continue
if gpu:
res = _round_one_gpu(sub, z_scale, tri, rz_push, sigma, close_r, open_r, dev)
else:
res = _round_one_cpu(sub, z_scale, tri, rz_push, sigma, close_r, open_r)
write = res & (out[slc] == 0)
out[slc][write] = lab
return out
# ---- 5a. shape-aware split [Yukimatsu.split_touching_shape_aware] ----
def _peak_markers(sub, dt, voxel_um, cfg):
fp = ellipsoid_se(*_radii_vox(cfg["min_peak_sep_um"], voxel_um, 1))
peaks = peak_local_max(dt, labels=sub.astype(np.uint8), footprint=fp, exclude_border=False)
if peaks.shape[0] > cfg["max_peaks"]:
peaks = peaks[np.argsort(dt[tuple(peaks.T)])[-cfg["max_peaks"]:]]
if peaks.size == 0:
return None, None
dt_max = float(dt[sub].max())
if dt_max <= 0:
return None, None
peaks = peaks[dt[tuple(peaks.T)] >= cfg["min_rel_dt_at_peak"] * dt_max]
if peaks.shape[0] < 2:
return None, None
markers = np.zeros(dt.shape, dtype=np.int32)
for i, p in enumerate(peaks, start=1):
markers[tuple(int(v) for v in p)] = i
return markers, dt_max
def _axes_ratio(mask):
props = regionprops(mask.astype(np.uint8))
if not props:
return 1.0
ev = np.sort(np.linalg.eigvalsh(props[0].inertia_tensor))
return float(np.sqrt(ev[-1] / ev[0])) if ev[0] > 0 else 1.0
def _merge_by_rules(parent, ws, dt, dt_max, cfg):
out = ws.copy()
total = int(parent.sum())
changed = True
while changed:
changed = False
ids = [i for i in np.unique(out) if i > 0]
if len(ids) < 2:
break
counts = {i: int((out == i).sum()) for i in ids}
dil = {i: binary_dilation(out == i) for i in ids}
to_merge = set()
for a in range(len(ids)):
for b in range(a + 1, len(ids)):
i, j = ids[a], ids[b]
touch = dil[i] & (out == j)
if not touch.any():
continue
if float(dt[touch].min()) >= cfg["neck_rel"] * dt_max:
to_merge.add((i, j)) # thick boundary -> one nucleus
continue
if any(counts[c] < cfg["min_child_voxels"] or
counts[c] < cfg["min_child_vol_ratio"] * total for c in (i, j)):
to_merge.add((i, j)) # tiny fragment
continue
if (_axes_ratio(out == i) < cfg["elong_merge_if_ratio_below"] and
_axes_ratio(out == j) < cfg["elong_merge_if_ratio_below"]):
to_merge.add((i, j)) # both children round -> false split
for i, j in to_merge:
if (out == i).any() and (out == j).any():
out[out == j] = i
changed = True
final = np.zeros_like(out)
for k, old in enumerate([i for i in np.unique(out) if i > 0], start=1):
final[out == old] = k
return final
def split_touching_shape_aware(masks, voxel_um=VOXEL_UM, cfg=None, progress=True):
"""Split labels whose distance transform has >=2 deep maxima separated by a thin neck.
The largest child keeps the parent label; others get new ids appended at the end."""
cfg = cfg or SPLIT_DEFAULT
out = np.asarray(masks).copy()
if out.max() == 0:
return out
bboxes = find_objects(out)
next_id = int(out.max())
it = range(1, len(bboxes) + 1)
it = tqdm(it, desc="shape-aware split", unit="obj", leave=False) if progress else it
for lab in it:
sl = bboxes[lab - 1]
if sl is None:
continue
sub = out[sl] == lab
if not sub.any():
continue
dt = distance_transform_edt(sub, sampling=voxel_um)
markers, dt_max = _peak_markers(sub, dt, voxel_um, cfg)
if markers is None:
continue
ws = watershed(-dt, markers=markers, mask=sub)
if ws.max() <= 1:
continue
ws = _merge_by_rules(sub, ws, dt, dt_max, cfg)
if ws.max() <= 1:
continue
counts = np.bincount(ws.ravel())
order = np.argsort(counts[1:])[::-1] + 1
out[sl][sub] = 0
out[sl][ws == order[0]] = lab
for p in order[1:]:
next_id += 1
out[sl][ws == p] = next_id
return out
# ---- 5b. cleanup [YukimatsuOnSynth2Stardist.cleanup_morphology] ----
def cleanup_labels(masks, voxel_um=VOXEL_UM, open_radius_um=CLEANUP_DEFAULT["open_radius_um"],
min_voxels=CLEANUP_DEFAULT["min_voxels"], progress=True):
"""Per label: largest connected component -> hole fill -> ellipsoid opening
(skipped when ``open_radius_um`` is 0/None) -> largest CC -> drop if < ``min_voxels``.
Returns (uint32 volume, stats dict)."""
masks = np.asarray(masks)
se = None
rz = ry = rx = 0
if open_radius_um:
rz, ry, rx = _radii_vox(open_radius_um, voxel_um)
se = ellipsoid_se(rz, ry, rx)
pad = max(rz, ry, rx) + 1
Z, Y, X = masks.shape
out = np.zeros(masks.shape, dtype=np.uint32)
kept = dropped = 0
bboxes = find_objects(masks)
it = enumerate(bboxes)
it = tqdm(it, total=len(bboxes), desc="cleanup", unit="obj", leave=False) if progress else it
for idx, sl in it:
if sl is None:
continue
lab = idx + 1
slc = (slice(max(0, sl[0].start - pad), min(Z, sl[0].stop + pad)),
slice(max(0, sl[1].start - pad), min(Y, sl[1].stop + pad)),
slice(max(0, sl[2].start - pad), min(X, sl[2].stop + pad)))
m = masks[slc] == lab
if not m.any():
continue
m = binary_fill_holes(largest_cc(m))
if se is not None:
m = largest_cc(binary_dilation(binary_erosion(m, structure=se), structure=se))
if m.sum() < min_voxels:
dropped += 1
continue
write = m & (out[slc] == 0)
out[slc][write] = lab
kept += 1
return out, dict(kept=kept, dropped=dropped, open_radius_um=open_radius_um,
min_voxels=min_voxels)
# ---- 5c. 2D -> 3D linking (hybrid) [HybridSegment.build_hybrid] ----
def link_2d_to_3d(masks2d, masks3d, overlap_thr=HYBRID_OVERLAP_THR, progress=True):
"""Hybrid volume: each independent 2D blob takes the label of the 3D object covering
the largest share of its pixels, if that share is >= ``overlap_thr``. Pixel boundaries
come from the 2D masks, z-linking from the 3D masks. Returns (uint32 volume, stats).
Not part of the current pipeline -- an alternative to 5d, claiming in the opposite
direction. Reachable via ``--hybrid``."""
if masks2d.shape != masks3d.shape:
raise ValueError("shape mismatch")
out = np.zeros(masks3d.shape, dtype=np.uint32)
st = dict(blobs=0, assigned=0, no_overlap=0, low_overlap=0)
Z = masks2d.shape[0]
it = tqdm(range(Z), desc="link 2d->3d", unit="slice", leave=False) if progress else range(Z)
for z in it:
s2, s3 = masks2d[z], masks3d[z]
for idx, sl in enumerate(find_objects(s2)):
if sl is None:
continue
m = s2[sl] == idx + 1
if not m.any():
continue
st["blobs"] += 1
under = s3[sl][m]
under = under[under > 0]
if under.size == 0:
st["no_overlap"] += 1
continue
counts = np.bincount(under)
best = int(counts.argmax())
if counts[best] / int(m.sum()) < overlap_thr:
st["low_overlap"] += 1
continue
out[z][sl][m] = best
st["assigned"] += 1
st["labels_out"] = int(np.unique(out).size - 1)
return out, st
# ---- 5d. 3D -> 2D anchoring (back-matching) [AnchorTo2D.py] ----
def anchor_to_2d(masks3d, masks2d, min_voxels=CLEANUP_ANCHOR["min_voxels"],
pad=ANCHOR_PAD_PX, progress=True):
"""Back-match a 3D label volume to the 2D masks, slice by slice.
For every 3D label and z-slice: if the label claims a 2D blob there (max-overlap claim,
conflicts between 3D labels resolved by overlap count), the slice becomes exactly that
2D blob; if the 2D segmentation saw nothing (caps, gap slices) the 3D slice is kept,
minus any 2D blobs claimed by other labels. Label ids are preserved.
Returns (uint32 volume, stats)."""
if masks2d.shape != masks3d.shape:
raise ValueError("shape mismatch")
objs = find_objects(masks3d)
labels = [i + 1 for i, sl in enumerate(objs) if sl is not None]
Z, Y, X = masks3d.shape
claims = {} # (z, label2d) -> (count, label3d)
it = tqdm(labels, desc="anchor: claims", unit="obj", leave=False) if progress else labels
for lab in it:
sl = objs[lab - 1]
for z in range(sl[0].start, sl[0].stop):
m3 = masks3d[z, sl[1], sl[2]] == lab
if not m3.any():
continue
under = masks2d[z, sl[1], sl[2]][m3]
under = under[under != 0]
if under.size == 0:
continue
l2, c = np.unique(under, return_counts=True)
for a, b in zip(l2, c):
key = (z, int(a))
if key not in claims or b > claims[key][0]:
claims[key] = (int(b), lab)
out = np.zeros(masks3d.shape, dtype=np.uint32)
st = dict(kept=0, dropped=0, slices_anchored=0, slices_kept3d=0, claims=len(claims))
it = tqdm(labels, desc="anchor: rebuild", unit="obj", leave=False) if progress else labels
for lab in it:
sl = objs[lab - 1]
y0, y1 = max(0, sl[1].start - pad), min(Y, sl[1].stop + pad)
x0, x1 = max(0, sl[2].start - pad), min(X, sl[2].stop + pad)
result = {}
vox = 0
for z in range(sl[0].start, sl[0].stop):
m3 = masks3d[z, y0:y1, x0:x1] == lab
if not m3.any():
continue
cp = masks2d[z, y0:y1, x0:x1]
under = np.unique(cp[m3])
under = under[under != 0]
mine = [int(l) for l in under if claims.get((z, int(l)), (0, -1))[1] == lab]
foreign = [int(l) for l in under if claims.get((z, int(l)), (0, -1))[1] != lab]
if mine:
res = np.isin(cp, mine)
st["slices_anchored"] += 1
else:
res = m3.copy()
if foreign:
res &= ~np.isin(cp, foreign)
st["slices_kept3d"] += 1
if res.any():
result[z] = res
vox += int(res.sum())
if vox < min_voxels:
st["dropped"] += 1
continue
for z, res in result.items():
tgt = out[z, y0:y1, x0:x1]
tgt[res & (tgt == 0)] = lab
st["kept"] += 1
return out, st
def anchor_and_finish(masks3d, masks2d, voxel_um=VOXEL_UM, cleanup=CLEANUP_ANCHOR,
rounding=ROUNDING_POST_ANCHOR, use_gpu=True, progress=True):
"""Anchoring followed by two finishing steps: an ellipsoid opening that removes shelves
one to three slices thick, and a rounding pass (``ROUNDING_POST_ANCHOR``) that closes the
residual slice steps."""
out, st = anchor_to_2d(masks3d, masks2d, min_voxels=cleanup["min_voxels"], progress=progress)
out, st_clean = cleanup_labels(out, voxel_um, progress=progress, **cleanup)
if rounding is not None:
out = round_caps(out, voxel_um, rounding, use_gpu=use_gpu, progress=progress)
st.update({f"cleanup_{k}": v for k, v in st_clean.items()})
return out.astype(np.uint32, copy=False), st
# ---- 5e. convenience chain ----
def refine(masks, voxel_um=VOXEL_UM, rounding=ROUNDING_DEFAULT, split=None,
cleanup=CLEANUP_DEFAULT, anchor_to=None, anchor_cleanup=CLEANUP_ANCHOR,
anchor_rounding=ROUNDING_POST_ANCHOR, use_gpu=True, progress=True):
"""compact labels -> round caps (``rounding``; None = skip) -> shape-aware split
(``split`` cfg; None = skip) -> cleanup -> optional anchoring to 2D masks
(``anchor_to`` = independent 2D label volume). Returns (uint32 volume, stats)."""
vol = relabel_compact_1n(masks)
stats = dict(labels_in=int(vol.max()))
if rounding is not None:
vol = round_caps(vol, voxel_um, rounding, use_gpu=use_gpu, progress=progress)
if split is not None:
vol = split_touching_shape_aware(vol, voxel_um, split, progress=progress)
stats["labels_after_split"] = int(np.unique(vol).size - 1)
vol, st = cleanup_labels(vol, voxel_um, progress=progress, **cleanup)
stats.update({f"cleanup_{k}": v for k, v in st.items()})
if anchor_to is not None:
vol, st = anchor_and_finish(vol, anchor_to, voxel_um, anchor_cleanup, anchor_rounding,
use_gpu=use_gpu, progress=progress)
stats.update({f"anchor_{k}": v for k, v in st.items()})
stats["labels_out"] = int(np.unique(vol).size - 1)
return vol.astype(np.uint32, copy=False), stats
# ============================================================================
# 6. SHAPE METRICS -- per-nucleus 3D and per-slice 2D
# [Compare3DvsMean2D.py metrics_3D / metrics_2D_sameMask / metrics_2D_indep;
# SliceSampling3Dvs2D.build_per_slice]
#
# Naming follows the committed per_nucleus.csv: ``R3D`` = Wadell sphericity psi
# (pi^(1/3) (6V)^(2/3) / A), ``S3D`` = min/max bounding-box extent;
# ``roundness_2d`` = 4 pi A / P^2 (circularity), ``sphericity_2d`` = min/max 2D bbox.
# ============================================================================
def _require_pandas():
if pd is None:
raise RuntimeError("the shape-metric tables need pandas; pip install pandas")
def metrics_3d(masks, voxel_um=VOXEL_UM, progress=True):
"""dict label -> {S3D, R3D, V_um3, A_um2, n_slices, z0, z1, y0, y1, x0, x1}."""
zsp, ysp, xsp = voxel_um
bboxes = find_objects(masks)
out = {}
it = enumerate(bboxes)
it = tqdm(it, total=len(bboxes), desc="metrics_3d", unit="obj", leave=False) if progress else it
for idx, sl in it:
if sl is None:
continue
lab = idx + 1
sub = masks[sl] == lab
if not sub.any():
continue
ext = ((sl[0].stop - sl[0].start) * zsp, (sl[1].stop - sl[1].start) * ysp,
(sl[2].stop - sl[2].start) * xsp)
V = float(sub.sum() * zsp * ysp * xsp)
try:
verts, faces, _, _ = marching_cubes(np.pad(sub.astype(np.uint8), 1), level=0.5,
spacing=voxel_um)
A = float(mesh_surface_area(verts, faces))
R3D = float(math.pi ** (1 / 3) * (6 * V) ** (2 / 3) / A) if A > 0 else float("nan")
except (RuntimeError, ValueError):
A, R3D = float("nan"), float("nan")
out[lab] = dict(S3D=float(min(ext) / max(ext)), R3D=R3D, V_um3=V, A_um2=A,
n_slices=int((sub.sum(axis=(1, 2)) > 0).sum()),
z0=sl[0].start, z1=sl[0].stop, y0=sl[1].start, y1=sl[1].stop,
x0=sl[2].start, x1=sl[2].stop)
return out
def blob_metrics_2d(blob):
"""(sphericity_2d, roundness_2d) of one boolean 2D blob; pixel units cancel."""
if not blob.any():
return float("nan"), float("nan")
yy, xx = np.where(blob)
by, bx = yy.max() - yy.min() + 1, xx.max() - xx.min() + 1
sph = float(min(by, bx) / max(by, bx))
props = regionprops(blob.astype(np.uint8))
if not props or props[0].perimeter == 0:
return sph, float("nan")
return sph, float(4.0 * math.pi * props[0].area / props[0].perimeter ** 2)
def metrics_2d_same_mask(masks, meta3d, progress=True):
"""Mean over the z cross-sections of the SAME 3D label (mode A)."""
out = {}
it = meta3d.items()
it = tqdm(it, desc="metrics_2d_same_mask", unit="obj", leave=False) if progress else it
for lab, m in it:
sphs, rnds = [], []
for z in range(m["z0"], m["z1"]):
blob = masks[z, m["y0"]:m["y1"], m["x0"]:m["x1"]] == lab
if not blob.any():
continue
s, r = blob_metrics_2d(blob)
if not math.isnan(s):
sphs.append(s)
if not math.isnan(r):
rnds.append(r)
out[lab] = dict(S_2Dmean_sameMask=float(np.mean(sphs)) if sphs else float("nan"),
R_2Dmean_sameMask=float(np.mean(rnds)) if rnds else float("nan"),
n_used_sameMask=len(sphs))
return out
def per_slice_bboxes(masks2d):
return [find_objects(masks2d[z]) for z in range(masks2d.shape[0])]
def metrics_2d_independent(masks3d, masks2d, meta3d, bboxes2d=None, progress=True):
"""Mean over the best-overlap INDEPENDENT 2D blob in each slice of the label (mode B)."""
bboxes2d = bboxes2d or per_slice_bboxes(masks2d)
out = {}
it = meta3d.items()
it = tqdm(it, desc="metrics_2d_independent", unit="obj", leave=False) if progress else it
for lab, m in it:
sphs, rnds, n = [], [], 0
for z in range(m["z0"], m["z1"]):
this = masks3d[z, m["y0"]:m["y1"], m["x0"]:m["x1"]] == lab
if not this.any():
continue
ids, counts = np.unique(masks2d[z, m["y0"]:m["y1"], m["x0"]:m["x1"]][this],
return_counts=True)
keep = ids != 0
ids, counts = ids[keep], counts[keep]
if ids.size == 0:
continue
best = int(ids[np.argmax(counts)])
if best - 1 >= len(bboxes2d[z]) or bboxes2d[z][best - 1] is None:
continue
blob = masks2d[z][bboxes2d[z][best - 1]] == best
if not blob.any():
continue
s, r = blob_metrics_2d(blob)
if not math.isnan(s):
sphs.append(s)
n += 1
if not math.isnan(r):
rnds.append(r)
out[lab] = dict(S_2Dmean_indep=float(np.mean(sphs)) if sphs else float("nan"),
R_2Dmean_indep=float(np.mean(rnds)) if rnds else float("nan"),
n_indep_matches=n)
return out
def per_nucleus_table(masks3d, masks2d=None, voxel_um=VOXEL_UM, progress=True):
"""DataFrame with one row per 3D label: 3D metrics + 2D-mean same-mask metrics and,
when ``masks2d`` is given, 2D-mean independent metrics."""
_require_pandas()
meta = metrics_3d(masks3d, voxel_um, progress)
same = metrics_2d_same_mask(masks3d, meta, progress)
indep = metrics_2d_independent(masks3d, masks2d, meta, progress=progress) \
if masks2d is not None else {}
rows = []
for lab in sorted(meta):
row = dict(label=lab, **{k: meta[lab][k] for k in ("n_slices", "S3D", "R3D", "V_um3", "A_um2")})
row.update(same[lab])
if indep:
row.update(indep[lab])
rows.append(row)
return pd.DataFrame(rows)
def per_slice_table(masks3d, voxel_um=VOXEL_UM, progress=True):
"""One row per (label, z) cross-section with the full 2D descriptor set: circularity,
bbox sphericity, ellipse aspect, eccentricity, solidity, convexity, Zheng-Hryciw
corrected roundness and centroid-radii sphericity (r_ins / r_circ)."""
_require_pandas()
xy = voxel_um[1]
bboxes = find_objects(masks3d)
rows = []
it = enumerate(bboxes)
it = tqdm(it, total=len(bboxes), desc="per_slice_table", unit="obj", leave=False) if progress else it
for idx, sl in it:
if sl is None:
continue
lab = idx + 1
for z in range(sl[0].start, sl[0].stop):
blob = masks3d[z, sl[1], sl[2]] == lab
if not blob.any():
continue
s, r = blob_metrics_2d(blob)
p = regionprops(blob.astype(np.uint8))[0]
area, perim = int(p.area), float(p.perimeter)
maj = float(p.axis_major_length)
aspect = float(p.axis_minor_length) / maj if maj > 0 else float("nan")
yy, xx = np.nonzero(blob)
dy, dx = int(yy.max() - yy.min()), int(xx.max() - xx.min())
try:
hull = ConvexHull(np.column_stack((yy, xx)))
convex = perim / hull.area if hull.area > 0 else float("nan") # 2D: hull.area = perimeter
except (QhullError, ValueError):
convex = float("nan")
if math.isnan(aspect) or math.isnan(r):
r_zh = float("nan")
else:
r_zh = r + (0.913 - (-0.528 * aspect * aspect + 1.259 * aspect + 0.270))
by_, bx_ = np.nonzero(blob & ~binary_erosion(blob))
d = np.hypot(by_ - yy.mean(), bx_ - xx.mean())
sph_r = float(d.min() / d.max()) if d.size and d.max() > 0 else float("nan")
rows.append(dict(label=lab, z_index=z, area_px=area, area_um2=area * xy * xy,
perimeter_um=perim * xy, roundness_2d=r, sphericity_2d=s,
aspect_2d=aspect, elong_2d=float(p.eccentricity),
eccbox_2d=dy / dx if dx > 0 else float("nan"),
convex_2d=convex, solid_2d=float(p.solidity),
roundness_zh_2d=r_zh, sphericity_radii_2d=sph_r))
return pd.DataFrame(rows).sort_values(["label", "z_index"]).reset_index(drop=True)
# ============================================================================
# 7. CLI
# ============================================================================
SYNTH_TAG = {"v1": "synth1", "v2": "synth2", "v3": "synth3"}
def _parse(argv):
p = argparse.ArgumentParser(
prog=os.path.basename(__file__), description=__doc__,
formatter_class=argparse.RawDescriptionHelpFormatter)
p.add_argument("--raw-dir", required=True,
help="folder of per-z tiffs (Zeiss ..._b0z<N>c0...tiff) or a single 3D tiff")
p.add_argument("--out-dir", required=True)
p.add_argument("--cellpose-weights", help="custom cpsam model; omit if M_indep.tiff exists")
p.add_argument("--cellpose-mode", choices=list(CELLPOSE_MODES), default="independent",
help="stage 1 front-end (default: independent = current pipeline)")
p.add_argument("--stardist-model-dir", help="StarDist basedir; omit if M_stardist_<v>.tiff exists")
p.add_argument("--stardist-python", help="python.exe of the stardist-gpu env; "
"omit to run stage 3 in this interpreter")
p.add_argument("--synth", choices=["v1", "v2", "v3"], default="v3")
p.add_argument("--variant", choices=list(VARIANTS), default="reduced")
p.add_argument("--split", action="store_true", help="also run the shape-aware split")
p.add_argument("--anchor", action="store_true", help="stage 5: back-match to the 2D masks")
p.add_argument("--hybrid", action="store_true",
help="additionally write the HybridSegment 2D->3D linking of the final masks")
p.add_argument("--metrics", action="store_true", help="stage 6: per_nucleus.csv + per_slice.csv")
p.add_argument("--mips", action="store_true",
help="write XY/XZ/YZ projections of the final masks (cast to uint16, so "
"label ids above 65535 wrap -- for visual QC only)")
p.add_argument("--voxel", default=",".join(str(v) for v in VOXEL_UM), help="Z,Y,X in um")
p.add_argument("--cpu", action="store_true", help="no GPU for cellpose and cap rounding")
p.add_argument("--no-progress", action="store_true")
p.add_argument("--force", action="store_true", help="recompute every stage")
return p.parse_args(argv)
def _stardist_stage(argv):
"""Stage 3 in its own interpreter: this file re-invoked with the stardist-gpu python."""
p = argparse.ArgumentParser(prog=f"{os.path.basename(__file__)} --stardist-stage")
p.add_argument("--stardist-stage", action="store_true")
p.add_argument("--input", required=True)
p.add_argument("--output", required=True)
p.add_argument("--model-dir", required=True)
p.add_argument("--model", default=STARDIST["model_name"])
p.add_argument("--tiles", default=",".join(str(t) for t in STARDIST["n_tiles"]))
p.add_argument("--prob", type=float, default=None)
p.add_argument("--nms", type=float, default=None)
a = p.parse_args(argv)
labels = run_stardist3d(imread_3d(a.input), a.model_dir, a.model,
tuple(int(t) for t in a.tiles.split(",")), a.prob, a.nms)
tifffile.imwrite(a.output, labels)
print(f"segmented {int(labels.max())} nuclei -> {a.output}")
return 0
def main(argv=None):
argv = sys.argv[1:] if argv is None else argv
if "--stardist-stage" in argv:
return _stardist_stage(argv)
a = _parse(argv)
voxel = tuple(float(v) for v in a.voxel.split(","))
if len(voxel) != 3:
raise SystemExit("--voxel needs three comma-separated values (Z,Y,X)")
progress = not a.no_progress
os.makedirs(a.out_dir, exist_ok=True)
P = lambda n: os.path.join(a.out_dir, n)
fresh = lambda path: a.force or not os.path.exists(path)
raw = load_raw(a.raw_dir)
print(f"raw {raw.shape} {raw.dtype}")
if not os.path.exists(P("raw_dapi.tiff")):
tifffile.imwrite(P("raw_dapi.tiff"), raw)
# --- 1. Cellpose 2D, independent -------------------------------------------------
if fresh(P("M_indep.tiff")):
if not a.cellpose_weights:
raise SystemExit("M_indep.tiff missing: pass --cellpose-weights")
kw = dict(voxel_um=voxel) if a.cellpose_mode == "3d" else {}
M_indep = CELLPOSE_MODES[a.cellpose_mode](raw, a.cellpose_weights, gpu=not a.cpu, **kw)
save_labels(P("M_indep.tiff"), M_indep)
else:
M_indep = tifffile.imread(P("M_indep.tiff"))
print(f"M_indep: {n_labels(M_indep)} 2D blobs ({a.cellpose_mode})")
# --- 2. synthetic DAPI -----------------------------------------------------------
synth_path = P(f"synth_dapi_{a.synth}.tiff")
if fresh(synth_path):
tifffile.imwrite(synth_path, build_synth_dapi(raw, M_indep, a.synth))
print(f"synth: {synth_path}")
# --- 3. StarDist 3D --------------------------------------------------------------
sd_path = P(f"M_stardist_{a.synth}.tiff")
if fresh(sd_path):
if not a.stardist_model_dir:
raise SystemExit(f"{os.path.basename(sd_path)} missing: pass --stardist-model-dir")
if a.stardist_python:
subprocess.check_call([a.stardist_python, os.path.abspath(__file__), "--stardist-stage",
"--input", synth_path, "--output", sd_path,
"--model-dir", a.stardist_model_dir,
"--model", STARDIST["model_name"]])
else:
labels = run_stardist3d(imread_3d(synth_path), a.stardist_model_dir)
tifffile.imwrite(sd_path, labels)
M_sd = tifffile.imread(sd_path)
print(f"stardist: {n_labels(M_sd)} labels")
# --- 4 + 5. Yukimatsu variant (+ optional split, + optional anchoring) -----------
rounding, cleanup = VARIANTS[a.variant]
tag = (f"M_{SYNTH_TAG[a.synth]}_yuki_{a.variant}"
+ ("_split" if a.split else "") + ("_anchored" if a.anchor else ""))
out_path = P(tag + ".tiff")
if fresh(out_path):
vol, stats = refine(M_sd, voxel, rounding=rounding,
split=SPLIT_DEFAULT if a.split else None,
cleanup=cleanup, anchor_to=M_indep if a.anchor else None,
use_gpu=not a.cpu, progress=progress)
save_labels(out_path, vol)
stats.update(synth=a.synth, variant=a.variant, split=a.split, anchor=a.anchor,
cellpose_mode=a.cellpose_mode, rounding=rounding, cleanup=cleanup,
voxel_um=voxel)
with open(P(tag + "_stats.json"), "w") as f:
json.dump(stats, f, indent=2, default=str)
print(json.dumps(stats, indent=2, default=str))
M_out = tifffile.imread(out_path)
print(f"final: {out_path} ({n_labels(M_out)} nuclei)")
if a.mips:
save_mips(P(tag), M_out)
# --- optional: HybridSegment 2D -> 3D linking of the final volume ----------------
if a.hybrid:
hyb_path = P(tag + "_hybrid.tiff")
if fresh(hyb_path):
hyb, st = link_2d_to_3d(M_indep, M_out, progress=progress)
save_labels(hyb_path, hyb)
print(f"hybrid: {hyb_path} {json.dumps(st)}")
# --- 6. shape metrics ------------------------------------------------------------
if a.metrics:
pn_path, ps_path = P("per_nucleus.csv"), P("per_slice.csv")
if fresh(pn_path):
per_nucleus_table(M_out, M_indep, voxel, progress).to_csv(pn_path, index=False)
print(f"metrics: {pn_path}")
if fresh(ps_path):
per_slice_table(M_out, voxel, progress).to_csv(ps_path, index=False)
print(f"metrics: {ps_path}")
return 0
if __name__ == "__main__":
sys.exit(main())

cbhs3d_onefile.py at commit 063985e, no license · at the source

Overview

Authors: Konstantin Hein^1,2,3,4,5,6, Humberto Romero‐Limon^1,2,3,6, Conrad Möckel^1,2,5, Anne Karasinsky⁴, Stephanie Möllmert^1,2, Jona Kayser^1,2,5, Alessio Zaccone⁷, Jochen Guck^1,2,5, Tomohisa Toda^1,3,4,6

ORCID iDs: Konstantin Hein, Humberto Romero‐Limon, Conrad Möckel, Anne Karasinsky, Stephanie Möllmert, Jona Kayser, Alessio Zaccone, Tomohisa Toda

Institutions: Friedrich-Alexander-Universität Erlangen-Nürnberg (Germany); Max Planck Institute for the Science of Light (Germany); German Center for Neurodegenerative Diseases (DZNE) Dresden (Germany); Universitätsklinikum Erlangen (Germany); Istituto Nazionale di Fisica Nucleare (Italy); University of Milan (Italy); Istituto Nazionale di Fisica Nucleare, Sezione di Milano (Italy)

Journal: Hippocampus, volume 36, issue 6, article e70139

Dates: received 2 July 2026; accepted 11 September 2026; published online 26 October 2026; in print November 2026

Type: Research article · Language: English

License: CC BY

Identifiers: DOI 10.1002/hipo.70139 · PMID 42934120 · PMCID PMC13710209 · OpenAlex W7223684705

Open access: hybrid, a free copy (OpenAlex)

Status: code verified (how the code is found and verified)

Categories: mouse (organism), cellular / molecular (subfield)

Methods: Connectivity, Statistics, Smoothing, state filtering, decompositions, Preprocessing, Machine learning, fMRI & imaging

Keywords: cellular packing order, hippocampus, quantitative image analysis, topology

MeSH: Hippocampus*, Image Processing, Computer-Assisted*, Neurons*, Animals, Astrocytes, Cell Nucleus, Male, Mice, Mice, Inbred C57BL, Software (* major topic)

Topic: Cell Image Analysis Techniques (Biophysics, Biochemistry, Genetics and Molecular Biology), according to OpenAlex

Funding: Deutsche Forschungsgemeinschaft (470322152 ‐ TO1347/3‐1, 497658532 ‐ TO1347/4‐1, 507965872 ‐ TO1347/5‐1, 563413271 ‐ TO1347/6‐1, 460333672 (CRC1540 EBM), 455449456); Schram Stifung; European Research Council (ERC‐2018‐STG, 804468 EAGER, ERC‐2023‐COG, 101125034 NEUTIME); Interdisciplinary Center for Clinical Research (IZKF) Erlangen (P162)

Citations: not cited yet (Europe PMC); 32 references in the paper

Abstract

The hippocampus is characterized by a stereotypical macroscopic structure, where the nuclei are densely and heterogeneously packed among different subregions of the hippocampus. Despite the fact that tissue‐specific cellular organization has been implicated in neural function, it has been technically challenging to quantitatively analyze mesoscopic cellular organization in the hippocampus due to its high cellular density. To overcome this technical hurdle, we developed Computational Biophysical Histomorphometry Software (CBHS), an automated image‐analysis pipeline aimed at quantifying nuclear shape and the order of the cellular ensemble in high‐density areas. When applied to the subfields of the hippocampus, we found that denser regions, most notably the dentate gyrus, were the most positionally, but least orientationally ordered. Nuclear shape exhibited a dependence on the local environment in a packing‐dependent manner. This association was cell‐type specific, with neurons, but not astrocytes displaying nuclear shape that varied with neighbor proximity, although astrocytes demonstrated greater intrinsic shape variance. The results reveal the presence of reproducible mesoscale cell packing order in hippocampal tissue, and are consistent with a nucleus‐driven mechanical coupling between neighboring cells. The present study provides a quantitative framework with which to understand mesoscopic tissue organization, thus enabling the formulation of testable hypotheses for future investigation.

Reproduced under the paper's license (CC BY), from the paper cited above.

Repository

Its files are read in the Code ↔ Paper reader above, with 3 matches between paragraphs and lines of code.

gitlab.mpcdf.mpg.de/konmi/image-analysis-ui

License: none: the authors keep all their rights

State: the link answers, verified on 1 November 2026

Evidence: files inventoried

Commit: 063985e1c2633ad5ccde4bfe29717844fb44e1e2, 18 September 2026

Languages: Python (7)

Size: 23 files, 7 scripts

Software Heritage: not archived

Found in: “Data Availability Statement”

Holds: README, environment (requirements.txt), continuous integration

Not found: license file, CITATION.cff, tests, documentation

Tools: NumPy (6 files), pandas (5 files), tifffile (4 files), scikit-image (3 files), SciPy (3 files), Cellpose (2 files), Matplotlib (2 files), CuPy (1 file), OpenCV (1 file), PyTorch (1 file)

Availability: 1 check, the latest on 1 November 2026: the link answers

1 November 2026: the link answers

8 files

code/2DAnalysis_Final_.py, Python, 543 lines
code/ExampleSHapes.py, Python, 158 lines
code/FunctionsClean.py, Python, 1,965 lines, 1 match
code/cbhs3d_onefile.py, Python, 1,250 lines, 2 matches
code/example_codes/ExampleAnalysis.py, Python, 99 lines
code/example_codes/ExampleSegmentation.py, Python, 40 lines
code/example_codes/ExampleTakeApart.py, Python, 67 lines
README.md, Text, 45 lines

At the source: gitlab.mpcdf.mpg.de/konmi/image-analysis-ui/-/tree/e1eea9740df5e9c3e5e0a185510e0f82145f9623/

The paper's code and data availability statement is in the Data section.

Tracing map

Proposed by the machine: these links were found in the paper and verified at the source, without human review. The map will receive a Zenodo DOI once one of the paper's authors has validated it with their ORCID.

What the map holds:

1 repository of the authors' code, each at its verified commit, with its license and how the link was found in the paper;
7 scripts, each with its path and the digest of its content;
3 matches between paragraphs of the paper and lines of the code (method lexical-v1);
neither the text of the paper nor the code itself.

Its JSON (tracing-map.json) is deposited on Zenodo with its DOI once the map is validated.

Data

No dataset and no data link were found in the paper.

Data Availability Statement

Images can be made available upon request to the leading author. All original code has been deposited at Gitlab under the permalink https://gitlab.mpcdf.mpg.de/konmi/image‐analysis‐ui/‐/tree/e1eea9740df5e9c3e5e0a185510e0f82145f9623/ (https://gitlab.mpcdf.mpg.de/konmi/image-analysis-ui/-/tree/e1eea9740df5e9c3e5e0a185510e0f82145f9623/) and is publicly available as of the date of publication. Any additional information required to reanalyze the data reported in this paper is available from the lead contact upon request.

Reproduced under the paper's license (CC BY), from the paper cited above.

Versions

The history of this record: each version stored by the harvester or made by a correction of its authors or of the maintainers of its code, and what changed in its facts. The texts of the paper (its abstract, its availability statements) are not part of it; versions that changed only those are not listed.

Version 1, 1 November 2026: the first record

Recorded: type, language, journal, volume, issue, pages, dates, 9 authors, 4 keywords, 10 MeSH terms, 4 funders, 28 references.

Cite

This paper

Hein, K., Romero‐Limon, H., Möckel, C., Karasinsky, A., Möllmert, S., Kayser, J., Zaccone, A., Guck, J., & Toda, T. (2026). Quantitative Assessment of Mesoscale Cellular Order and Organization in the Mouse Hippocampus. Hippocampus, 36(6), e70139. https://doi.org/10.1002/hipo.70139

BibTeX

@article{hein2026quantitative, author = {Hein, Konstantin and Romero‐Limon, Humberto and Möckel, Conrad and Karasinsky, Anne and Möllmert, Stephanie and Kayser, Jona and Zaccone, Alessio and Guck, Jochen and Toda, Tomohisa}, title = {{Quantitative Assessment of Mesoscale Cellular Order and Organization in the Mouse Hippocampus}}, journal = {Hippocampus}, year = {2026}, month = nov, volume = {36}, number = {6}, pages = {e70139}, publisher = {Wiley}, issn = {1050-9631}, doi = {10.1002/hipo.70139}, url = {https://doi.org/10.1002/hipo.70139}, pmid = {42934120}, pmcid = {PMC13710209} }

RIS

TY - JOUR AU - Hein, Konstantin AU - Romero‐Limon, Humberto AU - Möckel, Conrad AU - Karasinsky, Anne AU - Möllmert, Stephanie AU - Kayser, Jona AU - Zaccone, Alessio AU - Guck, Jochen AU - Toda, Tomohisa TI - Quantitative Assessment of Mesoscale Cellular Order and Organization in the Mouse Hippocampus T2 - Hippocampus J2 - Hippocampus PY - 2026 DA - 2026/11/01 VL - 36 IS - 6 SP - e70139 SN - 1050-9631 PB - Wiley DO - 10.1002/hipo.70139 UR - https://doi.org/10.1002/hipo.70139 LA - en ER -

CSL-JSON

{ "id": "10.1002/hipo.70139", "type": "article-journal", "title": "Quantitative Assessment of Mesoscale Cellular Order and Organization in the Mouse Hippocampus", "container-title": "Hippocampus", "author": [ { "family": "Hein", "given": "Konstantin" }, { "family": "Romero‐Limon", "given": "Humberto" }, { "family": "Möckel", "given": "Conrad" }, { "family": "Karasinsky", "given": "Anne" }, { "family": "Möllmert", "given": "Stephanie" }, { "family": "Kayser", "given": "Jona" }, { "family": "Zaccone", "given": "Alessio" }, { "family": "Guck", "given": "Jochen" }, { "family": "Toda", "given": "Tomohisa" } ], "container-title-short": "Hippocampus", "volume": "36", "issue": "6", "page": "e70139", "DOI": "10.1002/hipo.70139", "PMID": "42934120", "PMCID": "PMC13710209", "ISSN": "1050-9631", "publisher": "Wiley", "URL": "https://doi.org/10.1002/hipo.70139", "language": "en", "issued": { "date-parts": [ [ 2026, 11, 1 ] ] } }

The tracing map gets a citation of its own once an author has validated it and it has a DOI.

Similar papers

The papers with a page that share the most with this one: the tools found in their code, their categories, datasets, cited references and authors, the rarest counting most.

[1] doi:10.1038/s41586-025-09151-3 [code]: Morphodynamics of human early brain organoid development.
Journal: Nature
In common: Cellpose, CuPy, tifffile, 7 other tools, 1 reference
[2] doi:10.1113/ep093267 [code]: Spontaneous Ca<sup>2+</sup> signals in the developing mammalian cochlea of live mice under different anaesthetic regimes.
Journal: Experimental physiology
In common: Cellpose, CuPy, tifffile, 6 other tools, mouse, cellular / molecular, 1 reference
[3] doi:10.1038/s41467-024-55519-w [code]: In vivo spontaneous Ca<sup>2+</sup> activity in the pre-hearing mammalian cochlea.
Journal: Nature communications
In common: Cellpose, CuPy, tifffile, 6 other tools, mouse, cellular / molecular, 1 reference
[4] doi:10.1371/journal.pcbi.1014571 [code]: SynAPSeg: A novel dataset and image analysis framework for deep learning-based synapse detection and quantification.
Journal: PLoS computational biology
In common: Cellpose, tifffile, OpenCV, 6 other tools, 2 references
[5] doi:10.1038/s41587-022-01455-3 [code]: Scalable in situ single-cell profiling by electrophoretic capture of mRNA using EEL FISH.
Journal: Nature biotechnology
In common: Cellpose, tifffile, OpenCV, 6 other tools, mouse, cellular / molecular, 1 reference
[6] doi:10.1126/science.adq2084 [code]: Deep-tissue transcriptomics and subcellular imaging at high spatial resolution.
Journal: Science (New York, N.Y.)
In common: Cellpose, tifffile, OpenCV, 6 other tools, mouse, 1 reference
[7] doi:10.1016/j.cell.2022.12.010 [code]: Molecular and spatial signatures of mouse brain aging at single-cell resolution.
Journal: Cell
In common: CuPy, tifffile, OpenCV, 6 other tools, mouse, cellular / molecular, 1 reference
[8] doi:10.1016/j.isci.2024.109009 [code]: Spatial and single-nucleus transcriptomics decoding the molecular landscape and cellular organization of avian optic tectum.
Journal: iScience
In common: Cellpose, CuPy, OpenCV, 6 other tools, cellular / molecular
[9] doi:10.1073/pnas.2609132123 [code]: A human lysosomal storage disorder toolkit for decoding proteome landscapes in cortical-like and dopaminergic-like induced neurons.
Journal: Proceedings of the National Academy of Sciences of the United States of America
In common: Cellpose, tifffile, scikit-image, 5 other tools, cellular / molecular, 1 reference
[10] doi:10.1186/s40478-025-02048-2 [code]: Copper supplementation mitigates Parkinson-like wild-type SOD1 pathology and nigrostriatal degeneration in a novel mouse model.
Journal: Acta neuropathologica communications
In common: Cellpose, tifffile, OpenCV, 6 other tools, mouse

Contribute

The authors of this paper can claim it, correct its record and validate its tracing map, and the maintainers of its code (its owner, or a public member of its organization) correct what it says of their repository; anyone signed in can ask for its removal. Every request goes to OSCR's own machine, which answers it; your account page follows them.

Sign in with ORCID to claim this paper as one of its authors, correct its record or validate its tracing map: when the paper's metadata lists your ORCID iD, you are recognized at once. Maintainers of its code: sign in with GitHub, then claim the repository on your account page.

Add the badge to the code's README

If you are an author or a maintainer of this paper's code, you can add this badge. It says what OSCR knows of the paper (its code, its validated tracing map) and links the code to this page. Paste one of these snippets under the README's title. Keep the two marks around it, so that the badge can be replaced in place and is never added twice. OSCR changes nothing in your repository.

Markdown

<!-- oscr-badge:start -->
[![OSCR](https://openscicode.org/badge/doi/10.1002/hipo.70139.svg)](https://openscicode.org/badge/doi/10.1002/hipo.70139)
<!-- oscr-badge:end -->

reStructuredText

.. oscr-badge:start

.. image:: https://openscicode.org/badge/doi/10.1002/hipo.70139.svg
   :target: https://openscicode.org/badge/doi/10.1002/hipo.70139
   :alt: OSCR

.. oscr-badge:end

HTML

<!-- oscr-badge:start -->
<a href="https://openscicode.org/badge/doi/10.1002/hipo.70139"><img src="https://openscicode.org/badge/doi/10.1002/hipo.70139.svg" alt="OSCR"></a>
<!-- oscr-badge:end -->

If you choose to, you can have the badge kept up to date for you, with oscr badge --write from the command line or with the GitHub Action in your repository's workflows. Every other form and option (a card, AsciiDoc, LaTeX, shields.io) is on the badge page.

Request its removal

To ask OSCR to remove this record, the copies of its authors' scripts or its tracing map, use the removal request page: signed in, you say who you are, what to remove and why, then review and confirm the request. Published rules decide every request (how).

Discussion, reproductions, activity

Discussion: code errors and mismatches between the code and the paper about this paper and its code, from readers and its authors, beside the paper and its tracing map.

Report a mismatch between the code and the paper · Report a code error (signed in)

Reproductions: reports from readers who ran the authors' code and did not get the paper's result, what they ran, with which environment, commit and data.

Report a reproduction failure (signed in)

Activity: what happens around this paper: new versions of its record, its map's validation, discussions and reproductions. It opens with sign-in.