Search BANC EM with an LM colour-depth MIP (SRENs)

Goal

Take a 3-D light-microscopy image of a candidate SREN (sex-peptide- receptor rectal enteric neuron) and use it as a query against the BANC female adult connectome to recover the matching EM neuron(s). SRENs are a pair of cells at the posterior tip of the abdominal ganglion that express the sex peptide receptor and send efferent processes to the rectum; they are required for the rod-shaped fecal pellets (“RODs”) that mated females produce in response to seminal-fluid sex peptide. See Chae et al., Citations below, for the upstream/downstream circuit context.

The pipeline starts from the raw .lsm (a 40× confocal stack of the abdominal tip with the cord at an angle in the field of view) and lands a clean BANC ranking:

Extract GFP (label) and NC82 (anatomical reference) from the LSM at native voxdims via inst/python/lsm_to_nrrd.py.
PCA pre-rotate the LSM onto JRC2018VNCF axes (inst/python/lsm_pca_rotate.py). The NC82 tissue mask’s principal axes are aligned to template (X, Y, Z); the GFP CoG is anchored to a known abdominal-ganglion reference point so the rotated volume lands at the right tip in template coords. No Elastix optimiser run — the partial FOV / weak metric overlap made the rigid+affine optimiser unreliable, and PCA + an explicit anchor is more robust.
Mask GFP by the NC82 neuropil (inst/python/lsm_neuropil_mask.py). Over half of the raw GFP voxels (52% in our female sample) sit outside the cord neuropil — bleed-through into cuticle, rectal muscle, and background. Restricting to in-neuropil GFP triples the cMIP pixel-match count and promotes the right BANC matches.
Bridge the JRC2018VNCF NRRD into JRC2018VNCU_HR (the NeuronBridge VNC reference grid, 573 × 1119 × 219 voxels @ 0.461 × 0.461 × 0.7 µm) via the Saalfeld lab H5 displacement field shipped by navis-flybrains.
Render a colour-depth MIP of the bridged label channel with neuronbridger::nrrd_to_mip(target_space = "VNC").
Search that MIP against the BANC efferent MIP library (~1,033 flow == "efferent" cells in v888, 782 rendered) using colormip_search().
NBLAST the SREN intensity-weighted dotprops against the full efferent library using nat.nblast::nblast(..., normalised = TRUE) both directions and reporting the fwd+rev mean. Points are bridged from JRC2018VNCF µm into BANC µm via the full elastix JRC2018F→BANC chain (bancr::banc_to_JRC2018F(..., method = "elastix", inverse = TRUE)) — a stack of manual affine + elastix affine + coarse B-spline + fine B-spline. An earlier revision of this vignette used the simpler bancr::jrcvnc2018f_to_banc_tpsreg (a 5710-landmark TPS approximation) but that mapping collapses bilateral signal onto one side of the cord at the abdominal tip; the full elastix chain preserves it.
Render top-K BANC meshes overlaid on the SREN cloud as 2-D nat.ggplot panels in BANC native µm, plus a side-by-side cMIP of the SREN query and the top BANC hit, plus a bilateral-pair overlay for the current best hypothesis (EN00B016).

Why use BANC for this? BANC is the only female adult Drosophila CNS currently reconstructed at synapse resolution (v888 here, brain + VNC). Identifying the EM counterpart of the SRENs gives us their downstream connectivity — the descending input from SAG (Sex Peptide Abdominal Ganglion) neurons, the ascending afferents that signal back to the brain, and the rectal target. Colour-MIP search + NBLAST against BANC is the route to recovering these cells from a single LM volume.

Top-6 BANC efferent matches for the SREN query (cMIP, npmask)

What’s in the package

Neither the source .lsm (~480 MB, private) nor the BANC efferent MIP library (~12 MB across 782 PNGs) ship with neuronbridger. The vignette walks through the download steps; the only artefacts shipped are the rendered figures under inst/images/:

banc_colormip_sren_query.png — SREN query MIP in JRC2018VNCU_HR after PCA pre-rotation + neuropil mask.
banc_colormip_sren_top_hits.png — top-6 hits from colormip_search() (BANC EM mesh in purple, SREN LM cloud coloured by depth, faded by intensity). Each tile shows a top-down YX view + an XZ side projection. Plotted in BANC native µm; SREN cloud bridged via the full elastix JRC2018F→BANC chain.
banc_colormip_sren_nblast_top.png — same layout, ranked by NBLAST fwd+rev mean (normalised) against the full ~1,033-cell efferent library.
banc_colormip_sren_en00b016_pair.png — the two BANC cells whose manc_cell_type = EN00B016 (candidate SREN bilateral pair): 720575941680133053 (left, BANC cell_type = EFFabg09, ACh) + 720575941659397456 (right, BANC cell_type = EN00B016, octopamine).
banc_colormip_sren_vs_top_hit.png — side-by-side cMIP of the SREN query and the current top BANC cMIP hit (label written dynamically from live SeaTable annotations, since these have shifted over successive proofreading rounds).

A self-contained reproducer is at inst/scripts/banc_colormip_sren.R. The Python helpers it shells out to live under inst/python/: lsm_to_nrrd.py, lsm_build_nc82_mask.py, lsm_pca_rotate.py, lsm_neuropil_mask.py.

Prerequisites

# R packages
if (!require("remotes")) install.packages("remotes")
remotes::install_github("natverse/neuronbridger")
remotes::install_github("flyconnectome/bancr")    # BANC meshes + meta
remotes::install_github("natverse/nat.ggplot")    # 2-D neuron panels

# External tooling (one-off)
#   * Java 17+ (Saalfeldlab RenderTransformed JAR runtime)
#   * gsutil (Google Cloud SDK), authenticated to BANC's GCS bucket
#   * Python deps for the LSM pipeline:
#       conda activate r-reticulate
#       pip install tifffile SimpleITK scipy scikit-image pillow matplotlib
#
# The Saalfeldlab JRC VNC H5 transforms (~110 MB), via the navis
# `flybrains` Python package:
#
#   reticulate::py_run_string(
#     "import flybrains; flybrains.download_jrc_vnc_transforms()")
#
# This writes JRCVNC2018U_JRCVNC2018F.h5 to ~/flybrain-data/.

library(neuronbridger)
library(bancr)
library(nat); library(nat.ggplot); library(nat.nblast)
library(arrow); library(dplyr); library(ggplot2); library(patchwork)
library(scales)          # depth-coded SREN hue palette
library(reticulate)

Step 1 — From raw LSM to a JRC2018VNCF-aligned GFP NRRD

The raw LSM here is a partial FOV of the abdominal tip with the cord at an angle — the tip is roughly in the middle of the image and the rest of the VNC would be off to the top-right. A whole-VNC Elastix recipe doesn’t work (the optimiser keeps falling into local minima that put the cord on the wrong hemisphere), so we use a deterministic PCA pre-rotation instead.

WORK   <- tempfile("banc_colormip_work_"); dir.create(WORK, recursive = TRUE)
LSM    <- "/path/to/20250306_C24-53_EN1-5-1060_female_vnc_40x-2.lsm"
TPL    <- "/path/to/JRC2018_VNC_FEMALE_461.nrrd"
py     <- reticulate::conda_python("r-reticulate")
pydir  <- system.file("python", package = "neuronbridger")

# (a) Split the LSM into per-channel NRRDs. For this acquisition
#     ch0 = GFP, ch1 = unused red, ch2 = NC82.
system2(py, c(file.path(pydir, "lsm_to_nrrd.py"), LSM, WORK,
              "--prefix", "sren_female",
              "--gfp-channel", 0L, "--nc82-channel", 2L))

# (b) NC82 tissue mask (Otsu + close + largest CC + small dilation).
system2(py, c(file.path(pydir, "lsm_build_nc82_mask.py"),
              file.path(WORK, "sren_female_nc82.nrrd"),
              file.path(WORK, "sren_female_nc82_mask.nrrd")))

# (c) PCA pre-rotate onto JRC2018VNCF axes, anchor GFP CoG at
#     (131, 508, 92) um -- this is the GFP CoG of the upstream KDRC
#     registered TIF, taken as a known-good rough placement of the SREN
#     arborisation in template coords.
#
#     PC1 (longest mask axis) -> template Y (anterior-posterior)
#     PC2 (medium)             -> template X (lateral)
#     PC3 (shortest)           -> template Z (dorso-ventral)
#     Sign of PC1 is set so the GFP CoG projects positively onto PC1
#     (i.e. lands at high template Y, the posterior tip). det(R)=+1.
dir.create(file.path(WORK, "rot"), recursive = TRUE)
system2(py, c(file.path(pydir, "lsm_pca_rotate.py"),
              file.path(WORK, "sren_female_nc82.nrrd"),
              file.path(WORK, "sren_female_gfp.nrrd"),
              file.path(WORK, "sren_female_nc82_mask.nrrd"),
              TPL,
              file.path(WORK, "tip_mask.nrrd"),      # ignored, see --target-um
              file.path(WORK, "rot"),
              "--align-source", "gfp_cog",
              "--target-um",    "\"131,508,92\""))

The script writes rot/nc82_rot.nrrd, rot/gfp_rot.nrrd, and rot/nc82_mask_rot.nrrd, all sharing the same physical origin and voxdims (0.4 µm XY / 0.7 µm Z) with identity Direction. The R wrapper in inst/scripts/banc_colormip_sren.R then resamples these onto the JRC2018VNCF grid via SimpleITK (one-shot linear resample), so the output reads cleanly as a template-space NRRD.

Step 2 — Mask GFP by the NC82 neuropil

A simple Otsu of the resampled NC82 gives a neuropil mask; we zero any GFP voxel that falls outside it. About half of the raw GFP voxels in this acquisition sit outside the neuropil — fluorescence bleed-through into cuticle, rectal muscle, and background — so this step is a big lever:

system2(py, c(file.path(pydir, "lsm_neuropil_mask.py"),
              "nc82_PCAonly_JRC2018VNCF.nrrd",    # rotated + resampled NC82
              "gfp_PCAonly_JRC2018VNCF.nrrd",     # rotated + resampled GFP
              "gfp_PCAonly_neuropilmasked_JRC2018VNCF.nrrd",
              "--mask-out", "neuropil_mask_JRC2018VNCF.nrrd"))
#> NC82 Otsu threshold: 40.04
#> neuropil mask voxels: 5,462,433
#> GFP after neuropil mask: nonzero 36,798  (49.9% of pre-mask)

The masked GFP NRRD on the JRC2018VNCF grid is what feeds into the rest of the pipeline.

Step 3 — Bridge JRCVNC2018F → JRC2018VNCU_HR

ser_in_HR <- jrcvnc2018f_to_jrcvnc2018u_hr_h5(
  input  = "gfp_PCAonly_neuropilmasked_JRC2018VNCF.nrrd",
  output = file.path(WORK, "sren_female_in_JRC2018VNCU_HR.nrrd"))
# ~70 s on Apple silicon with 8 threads. Output: 573 x 1119 x 219
# voxels at 0.461 / 0.461 / 0.7 um.

The H5 (JRCVNC2018U_JRCVNC2018F.h5 from flybrains) stores its displacement field on the JRCVNC2018U grid, with dfield direction JRCVNC2018F → JRCVNC2018U (the file follows the brain <DEST>_<SRC>.h5 naming convention used by nat.jrcbrains). For F → U image resampling we want OUTPUT (U) → SOURCE (F) lookup — the inverse direction — so jrcvnc2018f_to_jrcvnc2018u_hr_h5() passes -i to RenderTransformed. Asking for the JRC2018VNCU_HR output grid directly works because JRC2018VNCU and JRC2018VNCU_HR share the same physical coordinate system; the JAR interpolates the dfield onto the requested grid in one shot.

Step 4 — Generate the query MIP

query_png <- nrrd_to_mip(
  input        = ser_in_HR,
  savefolder   = WORK,
  method       = "direct",
  target_space = "VNC",
  threshold    = 0.80,
  denoise      = "median3d",
  format       = "png",
  save         = TRUE,
  overwrite    = TRUE)

threshold = 0.80 keeps the top 20 % of non-zero voxels by intensity in the bridged volume; the neuropil mask already removed most of the dim halo so this floor is just a small extra noise cut.

SREN query colour-MIP in JRC2018VNCU_HR (neuropil-masked GFP)

Step 5 — Build the BANC efferent MIP library

LIB_DIR  <- file.path(WORK, "library", "mips")
dir.create(LIB_DIR, showWarnings = FALSE, recursive = TRUE)

# Shared root for every BANC GCS path used in this vignette.
GCS_ROOT <- "gs://lee-lab_brain-and-nerve-cord-fly-connectome"

# (a) BANC metadata — pull from the LIVE SeaTable `banc_meta` so we
# get the current `cell_type` + `manc_cell_type` bridge (both change
# with proofreading). The frozen `banc_888_meta.feather` on GCS is a
# snapshot; SeaTable is authoritative for annotation state and carries
# both the current `root_id` and the frozen `root_888` (which is what
# the cMIP filenames and precomputed layer object lists key on).
# Requires `BANCTABLE_TOKEN` in `.Renviron` — see `bancr::banctable_login`.
meta_live <- bancr::banctable_query(paste(
  "SELECT root_id, root_888, cell_type, manc_cell_type, cell_class,",
         "cell_sub_class, side, hemilineage, nerve, neuromere,",
         "neurotransmitter_predicted, peripheral_target_type,",
         "cell_function, flow",
  "FROM banc_meta"))
meta <- as.data.frame(meta_live) |>
  dplyr::filter(flow == "efferent") |>
  dplyr::select(root_888, side, hemilineage, nerve, neuromere,
                cell_class, cell_sub_class, cell_type, manc_cell_type,
                neurotransmitter_predicted, peripheral_target_type,
                cell_function)
nrow(meta)
#> [1] 1033

# (b) Match against the actual BANC color-MIP set (some efferents
# don't have a MIP rendered yet).
mips_root <- file.path(GCS_ROOT,
                       "neuron_colormips/template_alignment_240721",
                       "JRC2018_VNC_UNISEX_461") |> paste0("/")
all_mips  <- system2("gsutil", c("ls", mips_root), stdout = TRUE)
mip_ids   <- sub(".*/", "", sub("_in_JRC2018_VNC_UNISEX_461\\.png$", "", all_mips))
to_get    <- all_mips[mip_ids %in% meta$root_888]
length(to_get)
#> [1] 782

# (c) Parallel download (~2 MB total, takes <1 minute).
writeLines(to_get, file.path(WORK, "to_get.txt"))
system(sprintf("xargs -P 16 -I {} gsutil -q cp {} %s/ < %s",
               LIB_DIR, file.path(WORK, "to_get.txt")))

Step 6 — Search

res <- colormip_search(
  query        = query_png,
  library      = LIB_DIR,
  threshold    = 100L,    # channel-sum brightness cutoff (signal floor)
  z_tolerance  = 8L,      # depth-LUT-index tolerance
  xy_shift     = 3L,      # +/- 3 px translation grid
  mirror       = FALSE,   # SREN is unilateral
  mc.cores     = 8L,
  verbose      = TRUE)
res$root_888 <- sub(".*/(\\d+)_in_JRC2018.*", "\\1", res$path)

head(res[, c("root_888", "score", "n_match", "dx", "dy", "mirror")], 6)
#>             root_888     score n_match dx dy mirror
#> 1 720575941545785276 0.2122492     201  3 -3  FALSE
#> 2 720575941661415804 0.1774023     168  3 -3  FALSE
#> 3 720575941549868093 0.1636748     155  3 -3  FALSE
#> 4 720575941513645891 0.1499472     142  3 -3  FALSE
#> 5 720575941686473228 0.1488912     141  3 -3  FALSE
#> 6 720575941494055038 0.1478353     140  3 -3  FALSE

Neuropil-masking the query pushed the top cMIP score from ~0.072 (without the mask) to 0.21 here — three times more pixel matches because the noise floor that previously dominated the MIP is gone.

Step 7 — Cross-check with NBLAST fwd+rev mean (full elastix bridge)

We score the SREN dotprops against the full ~1,033-cell BANC efferent library using nat.nblast::nblast(..., normalised = TRUE) both directions and averaging the two — the standard NBLAST scorer. Points are bridged from JRC2018VNCF µm into BANC µm via bancr::banc_to_JRC2018F with method = "elastix" (the full manual-affine + elastix-affine + coarse-B-spline + fine-B-spline chain shipped at gs://lee-lab_brain-and-nerve-cord-fly-connectome/registrations/vnc_240721/), which preserves the bilateral SREN signal at the tip. Earlier revisions used bancr::jrcvnc2018f_to_banc_tpsreg (a TPS approximation of the same registration) and reported forward-only NBLAST because the TPS placed the bridged signal on one side of the cord — that shortcut is no longer needed.

library(nat.nblast)

# (a) Intensity-weighted SREN dotprops in BANC microns. Use ALL non-zero
# neuropil-masked voxels (no median filter, no largest-CC selection) so
# the bilateral pattern in the LSM survives.
ser <- nat::read.nrrd("gfp_PCAonly_neuropilmasked_JRC2018VNCF.nrrd")
storage.mode(ser) <- "integer"
idx <- which(ser > 0, arr.ind = TRUE)
ser_full <- data.frame(
  X = (idx[,1] - 1L) * 0.461122 + 0.461122/2,
  Y = (idx[,2] - 1L) * 0.461122 + 0.461122/2,
  Z = (idx[,3] - 1L) * 0.7      + 0.35,
  I = as.integer(ser[ser > 0]))

set.seed(11L); N <- 50000L
keep <- sample.int(nrow(ser_full), min(N, nrow(ser_full)),
                   prob = ser_full$I / sum(ser_full$I))
pts_F <- as.matrix(ser_full[keep, c("X","Y","Z")])
# Full elastix JRC2018F -> BANC (points in um).
pts_BANC_um <- bancr::banc_to_JRC2018F(
  pts_F, region = "vnc", banc.units = "um",
  inverse = TRUE, method = "elastix")
ser_dp <- nat::dotprops(pts_BANC_um, k = 5L)

# (b) Full-library NBLAST (~1,033 BANC efferents). The L2 dotprops are
# pre-cached locally as an .rds: fetch each cell's L2 skeleton with
# `bancr::banc_read_l2dp(meta$root_888)` once, drop NULLs/duplicates,
# and `saveRDS()` the resulting `neuronlist` so this step is offline.
banc_dps <- readRDS("banc_l2dp_efferent_full.rds")
banc_dps <- banc_dps[!sapply(banc_dps, is.null)]
banc_dps <- banc_dps[!duplicated(names(banc_dps))]

fwd <- nat.nblast::nblast(ser_dp, nat::as.neuronlist(banc_dps),
                          smat = nat.nblast::smat.fcwb, normalised = TRUE)
rev <- nat.nblast::nblast(nat::as.neuronlist(banc_dps), ser_dp,
                          smat = nat.nblast::smat.fcwb, normalised = TRUE)
nb <- data.frame(root_888    = names(fwd),
                 nblast_mean = (as.numeric(fwd) + as.numeric(rev)) / 2)
nb <- merge(nb, meta, by = "root_888")
nb <- nb[order(-nb$nblast_mean), ]
head(nb[, c("root_888", "nblast_mean", "cell_type", "manc_cell_type",
            "cell_sub_class", "side", "neurotransmitter_predicted")], 6)

Under this pipeline, the top-6 NBLAST hits are dominated by BANC cell_type = EFFabg40 and its close relatives (EFFabg43, EFFabg45), all sitting inside cell_sub_class = "abdomen_neurosecretory_cell" or "abd_Dsx_OA". These are doublesex-expressing sex-dimorphic abdominal neurosecretory cells, largely octopaminergic — a much better biological match for “sex-peptide-receptor rectal enteric neuron” than the motor-neuron interpretation an earlier revision of this vignette proposed. The MANC bridge (manc_cell_type) resolves these to legacy MANC names like MNad17, MNad48, MNad04, MNad69.

Top-6 BANC efferent hits ranked by NBLAST fwd+rev mean (full-elastix bridge)

Step 8 — Side-by-side: SREN query vs top BANC cMIP hit

SREN LM (left) vs top BANC cMIP hit (right)

Step 8b — The hypothesised SREN pair: manc_cell_type EN00B016

Neither the cMIP top nor the NBLAST top puts specifically the canonical SPR pair at rank 1, but both consistently rank a bilateral pair with manc_cell_type = "EN00B016" (exactly 2 cells, one per side in BANC) in the top ~10 % of the library. These are the biologically motivated SREN candidates:

left — root_888 720575941680133053, BANC cell_type = EFFabg09, MANC EN00B016, abdomen_motor_neuron, ACh (hemilineage 06A);
right — root_888 720575941659397456, BANC cell_type = EN00B016, MANC EN00B016, ventral_nerve_cord_neurosecretory_cell / abdomen_neurosecretory_cell, octopamine (hemilineage 00B).

They meet the bilateral-pair expectation for a sex-peptide-receptor population and land in the abd_Dsx_OA family the search prefers. Their MANC bridge to a single EN00B016 label suggests MANC groups them as one cell type across sides even though BANC currently annotates them as distinct cell_type names.

Hypothesised SREN pair (BANC manc_cell_type EN00B016)

What we recovered

Under PCA pre-rotation + neuropil mask + full-elastix JRC2018F→BANC bridge + NBLAST fwd+rev mean, the top-6 BANC efferents are all in the abd_Dsx_OA / abdomen-neurosecretory family (cell_class = ventral_nerve_cord_neurosecretory_cell), spread across a handful of BANC EFFabg* cell_types and their MANC bridges:

rank	BANC `cell_type`	MANC `manc_cell_type`	side	NT	cell_sub_class
1	EFFabg40	MNad17	right	octopamine	abdomen_neurosecretory
2	EFFabg40	MNad17	right	octopamine	abdomen_neurosecretory
3	EFFabg43	MNad48	right	octopamine	abd_Dsx_OA
4	EFFabg39	MNad04	left	octopamine	abd_Dsx_OA
5	EFFabg40	MNad17	right	dopamine	abdomen_neurosecretory
6	EFFabg40	MNad48	left	—	abdomen_neurosecretory

EFFabg40 (7 cells total in v888) accounts for most of the top hits; its MANC bridge is MNad17. The precise identity is most likely a member of the abd_Dsx_OA / abdomen-neurosecretory family; the best specific candidate for the SREN bilateral pair is manc_cell_type = EN00B016 (2 cells, one per side). Both L and R members appear in the top ~10 % of the ranking.

Step 9 — Cycling through candidates and inspecting the deformation

Two animated illustrations. Both are the kind of view Neuroglancer can show interactively but that a static vignette benefits from as a gif — one for the registration deformation itself, one for the candidate list.

SREN LM cloud with each top BANC candidate, in order

A slow cycle through the top-6 NBLAST hits + the two hypothesised SREN candidates (manc_cell_type EN00B016), overlaid on the elastix-aligned SREN LM cloud in BANC µm. This is the fastest way to “look through” the ranked list and see whether a given candidate’s dendritic tree traces the LM signal.

Top NBLAST candidates + hypothesised SREN pair on the SREN LM cloud (BANC um, full-elastix bridge)

The JRC2018VNCF -> BANC deformation, hull only, via `deformetricar`

The full-elastix chain (Step 7) is the recommended bridge; it does the heavy lifting from JRC2018F voxel intensity to BANC synapse density. As a separate experiment we tried an alternative registration using deformetricar — fitting a diffeomorphism between the two VNC neuropil hulls (the JRC2018VNCF outer surface and bancr::banc_vnc_neuropil.surf), driven by shape alone, no image content.

library(deformetricar)
# Source: JRC2018VNCF hull (isosurface of the female VNC template)
jrc_hull  <- Rvcg::vcgIsosurface(unclass(tpl >= thr), threshold = 0.5,
                                 spacing = c(0.461, 0.461, 0.7))
# Target: BANC VNC neuropil hull (bancr::banc_vnc_neuropil.surf, in um)
banc_hull <- bancr::banc_vnc_neuropil.surf

pre <- affine_prealign(jrc_hull, banc_hull, type = "similarity")
fit <- deformetrica_register(
  source = pre$apply(jrc_hull),
  target = banc_hull,
  kernel_width = 50,          # 50 um kernel
  object_type  = "SurfaceMesh",
  workdir      = "dfca_out",
  device       = "auto")

The flow gif shows the JRC hull morphing into the BANC hull; the top BANC NBLAST hit is drawn as the purple mesh for reference:

Deformetrica flow: JRC2018VNCF hull morphing onto the BANC VNC hull (top BANC hit in purple)

Bottom line: hull-only isn’t enough for this problem. The deformetrica warp fit hull → hull very well (~20 s on CPU), but when we apply the resulting deformation to the SREN cloud, the candidate ranks collapse (our hypothesised pair drops to rank ~700–800, and the top-15 is dominated by hind_leg_motor_neuron cells rather than abd_Dsx_OA neurosecretory ones). The reason is that the outer hull constrains only the perimeter; the interior gets whatever deformation the RKHS kernel prefers, which happens to shift the tip signal into the wrong internal compartment. To make deformetrica competitive here we’d need to feed it internal landmarks too — a family of matched neuropil sub-surfaces on both sides, similarly to what the mosquito-to-fly vignette does. The follow-up experiment below (MANC- cognate driven refinement) tries exactly that.

Following the pattern of the malecns-to-banc vignette in deformetricar, we tried to improve the JRC → BANC bridge in the abdominal region specifically by using cognate cell pairs from the male VNC connectome (MANC) that BANC’s SeaTable annotates via its manc_match column. There are ~102 abd_Dsx_OA / abdomen_neurosecretory_cell BANC cells with a MANC match.

Two headline pieces of work were needed to get this off the ground:

Bootstrap a MANC → JRCVNC2018F bridge. The nat.jrcbrains graph carries MANC ↔︎ JRCFIB2022M, JRCFIB2022M ↔︎ JRC2018M and JRCVNC2018F ↔︎ JRCVNC2018U/M, but no chain from MANC to any JRCVNC template. We built one from the cognate pairs themselves: for each of the 102 BANC-MANC pairs, bridge the BANC skeleton into JRC2018VNCF via the elastix chain and pair with the MANC skeleton in its native 8-nm-voxel space, then fit a thin-plate spline (nat::tpsreg) on the centroid landmarks.
Multi-object deformetrica fit with the malecns-to-banc pattern: deformetrica_register_multi() per (cell_type, side) using core-frame skeletons (Varifold), plus a weak whole-VNC-hull landmark anchor. Kernel = target-diagonal / 12 per the updated deformetricar README. Both sources and targets were placed in BANC µm (MANC bridged through JRC then through elastix) so the fit is a small refinement on top of the elastix bridge — that’s what the malecns pattern does for the brain.

library(deformetricar)

# ~102 BANC-MANC cognate pairs in abd_Dsx_OA / abdomen_neurosecretory
r <- bancr::banctable_query("SELECT root_id, cell_type, manc_match, side FROM banc_meta
                              WHERE (cell_sub_class='abd_Dsx_OA' OR
                                     cell_sub_class='abdomen_neurosecretory_cell')
                                AND manc_match IS NOT NULL AND manc_match != ''")

# Bootstrap a MANC -> JRCVNC2018F TPS from cognate-pair centroids.
# (BANC side goes native -> JRC via banc_to_JRC2018F(method="elastix");
#  MANC side comes from neuprintr in 8-nm voxels, scaled to um.)
# ... build src (MANC um centroid) + tgt (JRC um centroid) matrices ...
manc_to_jrc <- nat::tpsreg(sample = src_centroids, reference = tgt_centroids)

# Bring MANC skels into JRC then BANC (composed with the elastix bridge)
manc_in_banc <- nat::nlapply(manc_sk, function(x) {
  pts <- nat::xform(nat::xyzmatrix(x), manc_to_jrc)
  nat::xyzmatrix(x) <- bancr::banc_to_JRC2018F(
    pts, region = "vnc", banc.units = "um",
    inverse = TRUE, method = "elastix"); x
})

# Multi-object fit, per (cell_type, side) group, on core-frame skels
fit <- deformetrica_register_multi(
  srcs = manc_cores_by_group, tgts = banc_cores_by_group,
  kernel_width = diag / 12,          # ~50 um
  data_sigma = 0.5,
  landmarks = list(source = banc_hull_pts, target = banc_hull_pts),
  landmark_sigma = 5,
  timepoints = 10L, max_iterations = 40L, device = "cpu")

# Compose for the SREN cloud: elastix JRC -> BANC first, then refine.
pts_banc_elx <- bancr::banc_to_JRC2018F(sren_pts_jrc, region = "vnc",
                                         banc.units = "um",
                                         inverse = TRUE, method = "elastix")
pts_banc_refined <- nat::xform(pts_banc_elx, fit)

Results (scoring fwd + rev independently)

Under this refinement, the SREN’s hypothesised EN00B016 pair sits at rank_fwd = 67 / rank_rev = 57 for the right cell, and 90 / 84 for the left — essentially unchanged from the elastix-only baseline (58/59 for right, 88/86 for left). The refinement doesn’t materially move the candidate ranks.

But the top-10 under the refinement is much more biologically coherent than under elastix alone — nearly every top entry is a EFFabg40/EFFabg43/EFFabg45 cell in the abd_Dsx_OA or abdomen_neurosecretory_cell family, i.e. the correct broad identity for a sex-peptide-receptor abdominal neurosecretory population.

There is also a useful diagnostic that would be hidden by the fwd+rev mean: several top candidates by forward-NBLAST are far worse by reverse. E.g. the rank-1 forward hit (EN00B010) ranks 47 by reverse — a signature that the warped SREN cloud partially maps inside that neuron’s arbour but the neuron extends much further than the LM signal covers (or vice-versa). Scoring both directions separately catches this; the mean would flatten it. This is the reason the vignette now reports fwd and rev as independent columns.

Geodesic flow of MANC cognate neurons (pink) morphing onto their BANC counterparts (blue) inside the BANC VNC hull, under the deformetrica refinement fit

Verdict: The MANC-cognate refinement is a working malecns-to-banc-style pattern for the VNC and produces a top-hit list that’s more biologically consistent than the elastix bridge alone, but for this specific SREN identification it doesn’t rescue the candidates that the elastix bridge already puts in the top ~10 %. That’s an honest result: the elastix chain is already doing most of the alignment work for this cell family; the refinement is polish. It would matter more for cells with sparser cognate coverage or where the elastix bridge has bigger local errors — a fair use case for future runs.

Note also that a followup approach would be to keep going with the

Where the hull-only approach falls down (historical experiment)

The originally-experimented “just fit the two VNC outer hulls” case is documented above and remains as a cautionary tale — a nice demonstration of the READMEs “hull alone is not enough; add matched interior landmarks” rule. To skip through it: to make deformetrica competitive here we’d need to feed it internal landmarks too — a family of matched neuropil sub-surfaces on both sides, similarly to what the mosquito-to-fly vignette in deformetricar does. For now the full image-driven elastix chain is the better bridge; this experiment is a useful negative result and a template for a future improvement.

Step 10 — Open each match in Spelunker

Four Neuroglancer precomputed layers ship the SREN LM signal in BANC voxel space, corresponding to successive refinements of the JRC2018F→BANC bridge:

GCS path	registration	voxel size
`light_level/kdrc/SREN_female_aligned240721_to_BANC.ng`	bancr TPS (npmasked GFP)	800 nm
`light_level/kdrc/SREN_female_aligned240721_to_BANC_elastix.ng`	bancr TPS + per-specimen elastix	800 nm
`light_level/kdrc/SREN_female_aligned240721_to_BANC_elastix_400nm.ng`	bancr TPS + per-specimen elastix	400 nm
`light_level/kdrc/SREN_female_aligned240721_to_BANC_fullElastix_400nm.ng`	full elastix chain (recommended)	400 nm

Prefer ..._fullElastix_400nm.ng — that’s the registration-correct version, bilateral signal preserved, matching the NBLAST bridge used in Step 7.

LM_GS <- paste0("precomputed://", GCS_ROOT,
                "/light_level/kdrc/",
                "SREN_female_aligned240721_to_BANC_fullElastix_400nm.ng")
BASE  <- paste0("https://spelunker.cave-explorer.org/#!middleauth+",
                "https://global.daf-apis.com/nglstate/api/v1/",
                "6450802162925568")

top <- head(nb, 6)
top$ngl_url <- vapply(top$root_888, function(rid)
  bancr::banc_lm_scene(
    lm_url     = LM_GS,
    layer_name = "KDRC SREN (female, full elastix, 400 nm)",
    range      = c(1, 60),
    opacity    = 0.55,
    blend      = "additive",
    ids        = as.character(rid),
    url        = BASE,
    shorten    = TRUE),
  character(1))

Caveats

Anchor placement. The PCA rotation needs a target placement in template coords; we use the GFP CoG of an upstream KDRC-registered TIF (131, 508, 92 µm) as the anchor. If you don’t have that reference, use the fixed-mask CoG fallback in lsm_pca_rotate.py (default --align-source nc82_mask against a Y > 350 µm posterior tip mask) and accept a few µm of placement drift.
Sign disambiguation. The PCA orientation has a residual L/R ambiguity (PC2 sign). The default --pc2-sign +1 works for this acquisition; if your data lands on the wrong hemisphere, re-run with --pc2-sign -1.
Bridge accuracy — TPS is NOT enough at the abdominal tip. The first published revision of this pipeline used bancr::jrcvnc2018f_to_banc_tpsreg (a 5710-landmark TPS approximation of the elastix chain) to bridge JRC2018F→BANC. It measurably collapses bilateral signal to one side at the abdominal tip — the SREN X distribution in JRC2018VNCF was clearly bimodal around midline (X = 132 µm) but landed at X ∈ [308, 472] µm in BANC, entirely left of the BANC-tip midline (~X = 394 µm). We now use the full elastix chain via bancr::banc_to_JRC2018F(method="elastix"), which preserves the bilateral pattern. The bancr TPS is still fine for the brain bridge or for coarse-scale localisation, just not for the abdominal-tip topology.
JRC2018VNCU_HR bridge. The Saalfeldlab JRC VNC H5 (used for the cMIP rendering step) is generated from population averages and is small in the body of the VNC but can be imprecise at the periphery (where SREN axons leave). For sub-voxel precision, register your specific NC82 channel onto JRC2018VNCU_HR via Elastix multi-resolution affine + B-spline (the same recipe as lm_to_jrc2018u_elastix for the brain).
Library coverage. Only ~75 % of flow == "efferent" BANC v888 cells (782 / 1033) had a colour-MIP rendered at the time of writing. Re-run Step 5 as more get rendered.
NBLAST scoring. With bilateral signal preserved by the full elastix bridge, mean(fwd, rev) is the right scorer. The earlier forward-only recipe was a workaround for the TPS-collapsed unilateral query and is no longer needed.
Mirror. We still run cMIP and NBLAST with mirror = FALSE. The SREN GFP labelling is asymmetric (dominant left + dim right in this acquisition), so bilateral candidates surface as two individual unilateral hits (a left cell + a right cell of the same type family) rather than a single symmetric match. A mirror-aware second pass would lift a contralateral counterpart if you want the pairing explicit.
BANC ID drift and cell-type renames. Two annotation gotchas:
1. The BANC cMIP library and precomputed layers were rendered from the v888 snapshot (July 2024) and their filenames encode the frozen root_888. The live SeaTable’s default root_id column is the current post-proofreading ID and may differ. Always join by root_888 when going between the frozen assets (cMIP filenames, GCS banc_888_meta.feather) and the live SeaTable.
2. BANC cell-type names have shifted at least twice during this pipeline’s lifetime — a legacy “MNANmX” motor-neuron narrative gave way to MNadxxF, and the current annotation uses the EFFabg** prefix for VNC efferents with manc_cell_type as the bridge to MANC nomenclature. Prefer manc_cell_type as a reasonably stable pivot; treat any BANC cell_type as a moving target.
Cell-class re-assignment. The EFFabg40 / EFFabg43 / EFFabg45 cluster that dominates the top hits is annotated as cell_class = ventral_nerve_cord_neurosecretory_cell (specifically in the abd_Dsx_OA — doublesex-expressing sex-dimorphic — family) rather than as motor neurons. Earlier revisions of this vignette described these cells as “abdominal MNs”; that framing is superseded by the current curation.

Connectivity narrows the SREN candidate list

Independent of morphology, we can ask the BANC connectome directly: which efferent cells receive the strongest input, one or two hops away, from the curated sex-peptide sensory neurons (SPSNs)? If the SRENs express the sex peptide receptor and gate ROD production in response to mating, they should sit downstream of the SPSN pathway.

suppressMessages(library(bancr))
edges <- bancr::banc_edgelist(version = "v3", source = "gcs")
edges$pre <- as.character(edges$pre); edges$post <- as.character(edges$post)

spsn <- bancr::banctable_query(
  "SELECT root_id, side FROM banc_meta
   WHERE cell_sub_class = 'sex_peptide_sensory_neuron'")
# 6 curated SPSNs: reproductive-tract afferents, cell_class =
#   reproductive_tract_sensory_neuron, peripheral_target_type =
#   sex_peptide_receptor, cell_function = visceral_chemosensation.

# 1-hop: direct SPSN synapses.
downstream <- edges[edges$pre %in% spsn$root_id, ]
recip <- downstream |>
  dplyr::group_by(post) |>
  dplyr::summarise(spsn_in_syn = sum(count),
                   spsn_in_norm = sum(norm),
                   n_spsn_partners = dplyr::n()) |>
  dplyr::arrange(dplyr::desc(spsn_in_syn))

# 2-hop: SPSN -> intrinsic interneurons -> efferents, weighted by
#   path strength (norm_from_spsn * norm_from_intrinsic).
top30 <- head(recip$post, 30)
top30_wt <- setNames(head(recip$spsn_in_norm, 30), top30)
h2 <- edges[edges$pre %in% top30, ]
h2$path_wt <- h2$norm * top30_wt[h2$pre]

Two things pop out:

The top ~30 direct SPSN recipients are almost all intrinsic VNC interneurons (cell_class = ... unset or intrinsic; flow = intrinsic; cell_sub_class = ventral_nerve_cord_bilateral_ restricted / ipsilateral_restricted). Only one efferent abd_Dsx_OA cell sits in the top-15 (720575941518535347, EFFabg50 MNad12 left) — and it is explicitly annotated cell_function = "ovulation".
Ranking efferent abd_Dsx_OA / abdomen_neurosecretory_cell cells by 2-hop reach through those SPSN-driven intrinsics produces a very tight cluster whose top members are almost all annotated cell_function = "diuresis, stress_tolerance, pain_threshold". Diuresis matches Chae et al.’s rectal-ROD phenotype (RODs are the fluid-content-dependent fecal pellet morphology).

Top 10 connectivity + function-supported SREN candidates (BANC v888, by 2-hop SPSN reach, all cell_function = diuresis,stress_tolerance, pain_threshold):

rank	root_id	BANC cell_type	MANC	side	hemilineage	NT	2-hop weight
1	720575941575458056	EFFabg49	MNad13	left	00B	serotonin	0.0120
2	720575941582278265	EFFabg48	MNad03	left	06A	acetylcholine	0.0084
3	720575941652774421	EFFabg48	MNad03	left	06A	octopamine	0.0084
4	720575941550051873	EFFabg47	MNad52	left	13B	—	0.0083
5	720575941577258088	EFFabg44	MNad13	left	—	octopamine	0.0072
6	720575941594396267	EFFabg47	MNad03	right	13B	acetylcholine	0.0071
7	720575941482782659	EFFabg44	MNad13	right	—	octopamine	0.0071
8	720575941586627462	EFFabg47	MNad52	right	13B	gaba	0.0070
9	720575941563946038	EFFabg44	MNad13	left	—	acetylcholine	0.0070
10	720575941650150129	EFFabg49	MNad13	right	00B	gaba	0.0069

Bilateral pairs match the SREN morphology description. The EFFabg44, EFFabg47, EFFabg48, EFFabg49 cell types each contain matched left+right cells in this diuresis-annotated cluster.

This retracts the earlier EN00B016-pair hypothesis. The morphology + hull-registration analyses in the sections above put the hypothesised EN00B016 pair (720575941680133053 L + 720575941659397456 R) at rank ~59–90 by NBLAST fwd+rev — a plausible morphology score. But by connectivity:

The right cell (9397456) sits at abd-rank 18 by 2-hop SPSN reach (2-hop weight 0.006), with no functional annotation.
The left cell (0133053) has no 2-hop path from the SPSN cluster at all.

Neither cell is in the diuresis-annotated cluster. That’s not a tie-breaker in favour of EN00B016; it’s a signal against it, at least under the “SRENs are on the SPSN pathway and drive rectal diuresis” reading of Chae et al.

Recommendation for follow-up: focus on the EFFabg44 / EFFabg47 / EFFabg48 / EFFabg49 diuresis-annotated bilateral pairs as the current best SREN candidates. EFFabg49 (MNad13) is the top-scoring pair on 2-hop SPSN reach and worth inspecting first (left 720575941575458056, right 720575941650150129).

Caveat: the “diuresis, stress_tolerance, pain_threshold” cell_function label appears to be a family-level tag on the abd_Dsx_OA cluster; several cells outside the top-10 above (further down the 2-hop ranking) also carry it. Confirm-by-inspection is still needed — but the ranking + the functional label + the bilateral EFFabg47/48/49 pattern together are much stronger evidence than morphology alone.

Citations

SREN identification (this query): Chae HS, Subay OH, Kim DH, Yoon SE, Kim YJ. Rectal Enteric Neurons Optimize Fecal Production in Response to Mating Signals. Conference talk, APDNC4, 2026. Department of Life Sciences, Gwangju Institute of Science and Technology (GIST), Gwangju 61005, Republic of Korea (corresponding author: Young-Joon Kim, kimyj@gist.ac.kr; Korea Drosophila Resource Center, GIST). Identifies a pair of sex-peptide-receptor neurons at the abdominal-ganglion tip whose efferents innervate the rectum and gate ROD production in response to mating; downstream targets ascend to the brain via Trans-Tango.
BANC connectome: Bates et al. (2025). The Brain-And-Nerve-Cord fly Connectome (BANC): a complete adult female Drosophila central nervous system at synapse resolution. bioRxiv. https://www.biorxiv.org/content/10.1101/2025.07.31.667571v1
NeuronBridge / ColorMIP: Otsuna, Ito, Kawase. (2018). Color depth MIP mask search: a new tool to expedite Split-GAL4 creation. bioRxiv. https://doi.org/10.1101/318006
NBLAST: Costa, Manton, Ostrovsky, Prohaska, Jefferis. (2016). NBLAST: Rapid, Sensitive Comparison of Neuronal Structure and Construction of Neuron Family Databases. Neuron. https://doi.org/10.1016/j.neuron.2016.06.012
JRC2018 VNC templates: Bogovic et al. (2020). An unbiased template of the Drosophila brain and ventral nerve cord. PLOS ONE. https://doi.org/10.1371/journal.pone.0236495
flybrains Python package (Saalfeldlab JRC VNC H5 distribution): Schlegel et al., navis-org/navis-flybrains.