neu-mark¶

Annotations out of DVID, into columnar tables, and on into neuroglancer.

Two kinds of annotation:

point annotations — synapses, and anything else in a DVID annotation instance (labels_todo, bookmark_annotations).
body annotations — the per-body records in a keyvalue instance (labels_annotations): neuron name, status, who traced it.

It also writes a node’s ROIs as segmentation layers (below).

Tables come first and are the durable artifact: analysis reads them, and the neuroglancer layers are built from them, so a fetch is never repeated just to change how something is displayed.

Install¶

neuclease is a hard requirement for anything that talks to DVID and cannot be a pip dependency — it needs libdvid-cpp, vigra and dvidutils, all conda-only on flyem-forge:

mamba install -n neu-env -c flyem-forge -c conda-forge neuclease
pip install --no-deps -e ./neu-mark

Usage¶

Every fetch takes an explicit --bodies list, because the dataset has ~80M label ids and the great majority are single-voxel fragments nobody cares about. select-bodies is how you build that list.

Choosing the bodies¶

select-bodies ranks by synapse count, using the labelsz index DVID already maintains — AllSyn >= 10 returns 21,116 bodies in ~20 s on our dataset. Point --src at either the labelsz instance or the annotation instance it indexes (synapses), and its Base.Syncs is used to find the index.

Synapse count rather than voxel size, for three reasons: a body with no synapses is a fragment; many of the largest-by-volume bodies are glia; and scanning every label’s size is ~70 minutes (~80M labels at ~19k/s via fetch_listlabels, which pages serially by start = last_id + 1).

The threshold is on the total, not on pre and post separately. That is a domain constraint, not a simplification — sensory neurons may have no postsynapses, and a neuron projecting outside the traced volume may have no presynapses, so requiring both drops exactly the cells most worth looking at. --min-pre / --min-post exist when you do want that; both default to off.

Rank from a locked node. Body ids change with proofreading, so the list and the annotations should come from the same node — which is also why select-bodies writes a provenance record naming it, rather than a bare list of ids.

# what am I pointing at, and which node would I get?
neu-mark info --src dvid://dvid.example.org/93fdbc:main/synapses

# build the body list itself, from DVID's ranked synapse index (~20 s)
neu-mark select-bodies \
    --src dvid://dvid.example.org/93fdbc:main/synapses \
    --min-synapses 10 --out 'bodies_{uuid:8}.csv' --dvid-locked

# synapses for a body list -> points.parquet + relationships.parquet
neu-mark points \
    --src dvid://dvid.example.org/93fdbc:main/synapses \
    --bodies traced_neurons.csv \
    --out synapses_{uuid:8}/ --dvid-locked

# per-body annotations -> bodies.parquet  (a SEPARATE destination)
neu-mark bodies \
    --src dvid://dvid.example.org/93fdbc:main/labels_annotations \
    --bodies traced_neurons.csv \
    --out body_annotations_{uuid:8}/ --dvid-locked

--out is a directory for points and bodies, which write several tables, and the list’s file for select-bodies. Either may be a local path or s3://…, and may carry {uuid} / {uuid:N} / {uuid:full} / {branch} / {instance}, expanded from the resolved node. Naming an export after the node it came from matters because a branch ref means a different node tomorrow.

Tables¶

table	grain	key columns
`select-bodies --out` (the file named)	one row per selected body	`body`, `pre`, `post`, `syn`
`points`	one row per element	`body`, `z`, `y`, `x`, `kind`, `conf`, `user`, `roi`
`relationships`	one row per relationship	`rel`, `z/y/x`, `to_z/to_y/to_x`, `from_body`, `to_body`
`connections` (`--connections`)	one row per distinct (tbar, psd) pair	`pre_`, `post_`, `pre_body`, `post_body`
`bodies`	one row per body	`body`, `status`, `instance` (= neuron name), `user`, `json`

Coordinates are always named columns in z, y, x order. tables.positions_zyx and tables.positions_xyz are the only places a positional coordinate array is built, and therefore the only places the order is decided — a mirrored annotation is a valid annotation in the wrong place, and nothing else would catch it.

parquet is the default format because csv cannot preserve an integer column that contains nulls, and to_body is one by design — an unresolved partner is null. csv has no types, so pandas infers on read, and a column mixing integers with blanks never infers as an integer: measured here it returns as str, on older pandas as float64, which silently rounds ids above 2^53. A non-nullable uint64 column does survive csv intact, which makes this the worse failure mode — the corruption is confined to the one column nobody thinks to check. --format csv is available and warns, naming the columns that will not survive.

Which neuropil each synapse is in¶

--rois (or points(..., rois=[...])) adds a roi column to the points table, via neuclease’s determine_point_rois. It samples a combined ROI volume at scale 5, so it is cheap, and it takes x/y/z columns — which our tables have by name, so there is no flip and no chance of one.

neu-mark points --src @synapses --bodies traced.csv --out syn/ \
    --rois 'ME(L),ME(R),LO(L),LO(R),AGNG,PGNG,SNP(L),SNP(R)'     # or --rois @neuropils

--rois SPEC.toml:LAYER takes the ROIs of a roi-layers layer, so synapses are labelled with the same set the viewer shows: rois.toml:neuropil-parts for the parts of a neuropil tree.

The ROIs are read from --src’s node unless --roi-source dvid://SERVER/UUID names another, for ROIs redrawn on a later node: the points and their body ids still come from --src, so they keep matching a segmentation exported from it. The provenance record names the ROI node.

Only the points table gets it. A relationship spans two points that may be in different neuropils, so one column on it would have to pick a side; a join back to points answers either. tables.body_roi_counts(points) aggregates per body and kind.

For segment-property rules, segment-properties --synapses syn/ gives each body a synapse_neuropils column: the share of its pre, post and all synapses in each ROI, {"pre": {"ME(L)": 0.62, ...}, "post": {...}, "all": {...}}, largest first. With --roi-spec rois.toml, synapse_neuropils_top groups the same shares by the spec’s top-level neuropils. A rule can then tag a cell by where it connects rather than by its name:

@rule(needs=("synapse_neuropils_top",))
def region(r):
    """where most of the cell's synapses are"""
    shares = r["synapse_neuropils_top"]["all"]
    return next(iter(shares), None) if shares and max(shares.values()) > 0.5 else None

The ROI set is required — there is no “all”. The combined volume is built by writing each ROI in turn, so passing every ROI on the node would attribute a point in ME(L) to whichever of ME(L) / OL(L) / all_neuropils was written last. A name that is not an ROI instance on the node is caught before any fetch, with close matches suggested.

Where ROIs do intersect the later one wins, and there is no principled tie-break — so this proceeds and measures it, in the unit that decides whether to care: the volume is unpacked a second time with the order reversed (the ROI fetch is not repeated, so this is nearly free) and the two labellings compared. The report says how many of your synapses change hands:

rois: 989 of 1919 synapses inside one of 3 ROIs (51.5%), 930 outside every one
  OL(L) (919), AGNG (70)
  1 ROI pairs intersect; 919 synapse(s) (47.89%) sit in an intersection and are
  attributed by ROI ORDER alone — {'OL(L) | ME(L)': 919}

--strict-rois refuses instead, for a strict partition; the dataset has deliberately subtracted variants (INP(-ATL)(L), PENP(-AMMC), VLNP(-AOTU)(L)) for that.

The partner-resolution caveat¶

A DVID relationship points at a coordinate, not a body. The partner’s body is recovered by joining that coordinate against the elements table — which works because elements are fetched per body, so every row already knows its own body. But an edge only resolves when both endpoints were fetched, so the yield depends on how much of the connectome --bodies covers. Measured on our dataset, top-N bodies by presynapse count:

bodies	points	distinct pairs	both ends resolved
20	6,148	23,652	0.2%
100	26,748	75,901	2.0%
400	80,479	226,334	4.6%

So the commands report the match rate as a headline number, and unresolved partners are kept as nulls rather than dropped. --drop-unmatched exists for once the rate is high: dropping earlier turns “the body list did not cover this partner” into “this synapse does not exist”, which at low coverage presents a mostly-incomplete connectome as a complete one.

Every synapse: `points --all`¶

neu-mark points --src @synapses --all --out syn/ --connections --dvid-locked

--all removes the caveat. It fetches the whole instance by region (DVID’s /blocks, through neuclease’s fetch_synapses_in_batches) and looks up every element’s body in the synced labelmap at the same node (fetch_bodies_for_many_points: the supervoxel under each point, then DVID’s whole supervoxel-to-body map), so both ends of every pair are known. Measured on a 2.3M-element instance: the fetch takes about 30 s, and 6.6% of its pairs touch no body of a ≥10-synapse list, so a body list misses them entirely. An element on background label has body 0, which the connections table treats as no body. The points table keeps each element’s sv, the supervoxel, which makes a later lookup cheaper.

The tables are the same shape as for a body list. annotation-source --listed BODIES then marks which ends are on the bodies the segmentation shows (pre_listed / post_listed), so the viewer can hide the synapses between fragments.

Reading the whole instance¶

neu-mark bodies --src @bodies --all --out 'all_bodies_{uuid:8}/' --dvid-locked

This never calls /keys, which is unreliable at this size — 58,394 keys took 52 s twice and then returned 504 Gateway Time-out from the proxy in front of DVID, and there is no count endpoint to use instead. It reads a cover of the key space as bounded keyrangevalues requests, so completeness comes from the cover being exhaustive rather than from knowing a total. Values cost essentially nothing over keys (56.7 s against 54.3 s).

Two properties of DVID worth knowing if you use the underlying functions directly: its key ranges are inclusive at both ends, so consecutive ranges overlap by one key and summing per-range counts overcounts; and a range that is too big to finish inside the proxy’s window is subdivided rather than retried, because repeating it unchanged just fails again more slowly.

In a notebook¶

The same fetches, returning DataFrames and writing nothing:

from neu_mark import source, select_bodies, points, body_annotations

src = source("@synapses", locked=True)        # node pinned once
sel = select_bodies("@counts", min_synapses=10, locked=True)   # body, pre, post, syn
pts, rels = points(src, sel.head(50))          # bodies: frame, Series, list or a path
ann = body_annotations("@bodies", sel.head(50), locked=True)

neu_mark.ops is the other half — same fetches, but writing tables and a provenance record to a destination. Both go through the same functions, so what you see in a notebook is what a run would write.

Site defaults (optional)¶

An uncommitted TOML file supplies the server, uuid and instance names, so @synapses stands in for a full URL. It is a URL builder, not a fallback: nothing consults it implicitly, @name is required so a saved command visibly depends on it, and both the CLI and the library print what a reference resolved to.

It is found as neu-mark.toml or configs/neu-mark.toml in the working directory or any parent, so one file at the workspace root serves every repo and every notebook subdirectory beneath it — the way pyproject.toml is found. $NEU_MARK_CONFIG overrides. There is deliberately no machine-wide location: a hidden ~/.config file applies to every shell on the host and is the kind of invisible state that makes one command behave differently for two people.

# neu-mark.toml, at the workspace root
[dvid]
server = "dvid.example.org"
uuid   = "93fdbc:main"
locked = true

[instances]
synapses = "synapses"
bodies   = "labels_annotations"
counts   = "synapses_labelsz"

[roi_sets]
neuropils = ["ME(L)", "ME(R)", "LO(L)", "LO(R)", "AGNG", "PGNG"]

$ neu-mark select-bodies --src @counts --min-synapses 400 --out 'sel/{uuid:8}.csv' --dvid-locked
--src @counts  ->  dvid://dvid.example.org/93fdbc:main/synapses_labelsz
--out sel/{uuid:8}.csv  ->  sel/821d68d2.csv

Inspecting the annotation strings¶

instance is the field that carries the information and it is dirty in bounded ways. neu_mark.explore is a set of functions for looking at it — DataFrame in, DataFrame out, no printing and no I/O, so it is meant for a notebook:

from neu_mark import explore as ex

ex.instances(bodies, "instance")                     # distinct strings -> body counts
ex.tokens(bodies, "instance", drop=r"\((L|R)\)")     # the suffix vocabulary, side removed
ex.variants(bodies, "instance")                      # what normalization repaired ('truncated?')
ex.near(ex.tokens(bodies, "instance"), "fragment")   # fuzzy candidates for misspellings

near and variants catch different things and you want both: whitespace and a trailing ? are repaired by normalize, so truncated? never appears in the token histogram and near finds nothing left to fix. variants is what shows you the repair happened.

Once you have a guess about the vocabulary, try it out. Here a rule is nothing more than a function of one row, and rules is a plain dict of them:

guesses = {"cell_type": lambda r: (r["instance"] or "").split("_")[0] or None}

ex.apply_rules(bodies, guesses)                  # one column per rule
ex.coverage(bodies, guesses)                     # bodies, coverage, distinct, top
ex.unparsed(bodies, guesses, "cell_type", column="instance")   # what it missed, by bodies
ex.compare(bodies, "type", guesses["cell_type"], column="instance")   # curated vs your parse

A rule takes one row (a Series, with missing values as None) and returns a scalar, a sequence for genuinely multi-valued facets — column labels are, since some central complex neurons innervate several — or None for “did not fire”, which stays visible rather than becoming an error.

Writing rules¶

A dict of lambdas is right for the first five minutes and wrong for anything you publish from. neu_mark.rules is the same contract with the declarations added, and a dataset-specific RuleSet is loaded by path from a module living with the dataset. Nothing is built in. The segment-properties tutorial walks through a whole module.

# wasp_rules.py
import re
from neu_mark import explore as ex, rule

LABEL = "instance"                   # the field naming each body; rules read it by default
KEEP = ["instance", "type"]          # explicit; nothing else reaches a viewer

@rule
def side(r):
    """Hemisphere, from the parenthesized (L)/(R)."""
    m = re.search(r"\((L|R)\)", str(r["instance"] or ""))
    return m.group(1) if m else None

@rule(multi=True)
def column(r):
    """Column labels; genuinely multi-valued."""
    return re.findall(r"_([A-Z]\d+)(?=$|[_(])", str(r["instance"] or ""))

@rule(drop=True)
def noise(r):
    """Bodies this dataset marks as not real data."""
    tokens = {t.lower() for t in ex.split_tokens(r["instance"])}
    return next((t for t in ("irrelevant", "block") if t in tokens), None)

from neu_mark import rules
RS = rules.from_module("wasp_rules.py")

What a declaration buys you, beyond the bare function:


`needs`	the columns it reads, checked up front — a missing one is an error naming the rule, not an `AttributeError` 8,000 rows in. Left out, the rule reads the module’s `LABEL` field
`multi`	checked, not decorative: a single-valued rule returning two values raises, which catches `re.findall` where `re.search` was meant
`prefix`	the facet namespace, defaulting to the rule’s name — `glomerulus` gives `glomerulus:VA1v` with nothing declared. `prefix=""` mints a bare tag, right for a flag. The colon is added for you
`drop`	the rule excludes the body it fires on and returns the reason, which is counted in the report
`consumes`	which tokens of the instance this rule accounts for, so `RuleSet.remainder` can hand you what no rule recognized. Derived from the returned values by default — declare it only where a rule normalizes an alias (`NCL` → `nucleated`)
`tag=False`	produces something other than a tag. `label` and `description` are reserved and set those properties; anything else is computed, reported and written as a column of the `--table` output, but not put in the document
the docstring	becomes the description — in the coverage report, and on every tag the rule mints, so the viewer explains itself

Testing one rule, or all of them, on one string:

RS["side"].test("Tm2_A2(L)")        # 'L'
RS["side"].test("AGNG")             # None — did not fire
RS.explain("LN_C5(L)_NCL")          # {'side': 'L', 'column': ('C5',), ...}

print(RS)                           # each rule's repr: declarations, then what it is for
RS.describe()                       # the same as a frame, to sort or filter
RS.coverage(bodies)                 # what they did to your data, plus what they declare
RS.unparsed(bodies, "column")       # what to fix next, ranked by bodies

test and explain build a synthetic row through the real path, so missing fields arrive as None rather than pd.NA — the difference that makes r["instance"] or "" work.

What no rule recognized¶

RS.remainder("LN_C5(L)_NCL")     # 'LN'
RS.consumed("LN_C5(L)_NCL")      # {'side': ('L',), 'nucleated': ('NCL',), 'column': ('C5',)}

remainder strips every token the rules accounted for and rejoins the rest, which collapses runs of separators and strips them from the margins. It is the raw material for a cell-type parse that does not carry its own copy of every vocabulary the rules already encode — the failure mode being that you add a flag to a rule, forget to add it to the parser, and the flag starts appearing inside cell-type names.

Removal is token-granular rather than substring-based, because side returns "L" and deleting that substring from LC10_C5(L) would eat the L of LC10. Rules whose value is the matched token need no consumes; only the alias-normalizing ones do.

Be aware what comes back is a residue, not a cell type: on this dataset roughly three quarters of non-empty remainders are a bare neuropil name (AGNG, VLNP) — a body nobody has identified — so making a facet out of one needs a vocabulary that can tell those apart.

Publishing rules as segment properties¶

segment-properties runs the rules of the module --rules names: its RULES, or every @rule at module level in declaration order. Nothing is built in — not side, not a fragment flag — because which part of a name means what is a fact about your dataset, and another dataset’s vocabulary would mint plausible tags that fire on almost nothing. Without --rules the document carries labels and numbers only. The label is the module’s LABEL field (or --label FIELD); with neither, segments are known by id. --drop RULE excludes the bodies a rule fires on instead of tagging them, using the module’s own definition.

$ neu-mark segment-properties --src @bodies --dst ./seg --bodies sel/bodies.csv \
      --rules wasp_rules.py --dry-run
segments: 21116   tags: …   of 2111 allowed   document: … KiB
rules: wasp_rules.py
excluded: 394 — block (7), chunk (4), irrelevant (380), unknown (3)
  noise             394 bodies (  1.8%)      4 distinct   [drop]
  side            16238 bodies ( 76.9%)      2 distinct
  fragment         7109 bodies ( 33.8%)      1 distinct
  truncated        6625 bodies ( 31.5%)      1 distinct
  column               …                                          <- yours
DRY RUN — nothing written. ./seg is untouched.

(Counts are from the ≥10-synapse body list; the elided ones depend on your rules.)

--meshed-bodies in place of --bodies takes the bodies that have a mesh in --dst, so the properties cover what neu-morpho run meshed into that volume. It reads unsharded meshes only; for a sharded source, pass --bodies.

The tag vocabulary is written in rule order, because the viewer lists tags in vocabulary order: each facet’s values together, facets in the order their rules are declared, bare flags last, and within a facet case is ignored and numbers compare as numbers (col:A2 before col:A10). In the Segments tab’s query box, #side:L keeps a tag and -#fragment drops one, /Tm[0-9] matches the label as a regular expression, syn>=100 constrains a number, >syn sorts by it descending (< ascending), and |post adds it as a column.

Run --dry-run first: it builds and reports without writing anything, which --no-link does not — that still writes the document. A rule that fires on 2% of bodies looks exactly like a facet the data does not have, and only the coverage number tells them apart.

--table FILE writes the same build as a table, a row per requested body, dry run or not:

neu-mark segment-properties --src … --dst … --bodies bodies.csv --rules rules.py \
    --dry-run --table work/segprops.parquet

t = pd.read_parquet("work/segprops.parquet")
t[t["tag:type"].isna() & t.excluded.isna()][["body", "instance"]]   # what no rule typed

Each row has label, the full tags list, a tag:<prefix> column per facet (a list for a multi rule, else a value or null), a flag:<tag> boolean per bare tag, the counts, and then the record the rules read (instance, type, …; dict values such as the synapse neuropil shares as JSON text). Bodies left out of the document stay in the table, with the reason in excluded: the drop rule’s, or no record. Parquet carries the run’s provenance in its own metadata (neu_mark.io.read_embedded_provenance), and a <name>.provenance.json lands beside it.

Two guards, because a bad rule produces a valid source that no viewer can use.

--max-tags caps the vocabulary and takes either form: below 1 it is a fraction of the segments written, 1 or above an absolute count. The default is 0.10, because the honest size is dataset-dependent — 11,752 distinct type values over ~165k male-CNS bodies is 7% and entirely reasonable, while the same number over 21k bodies means a rule is handing back the instance string. The error names the rule with the most distinct values.

The branch is chosen by value, never by whether you wrote a decimal point: 1 and 1.0 are the same number to TOML and to JSON, so a setting meaning “all of them” would come back from a config meaning “one tag”. A converted fraction is floored at 100 tags so a short list is not judged by a ratio; an explicit count is a statement and is respected as given, floor or no floor.

Distinct facets are capped separately, which catches the other failure: tags that are each their own field can be a small vocabulary and still leave nothing to group on.

Because the rules now decide what the source says, the provenance sidecar records them — the resolved path, a sha256, git details if the module happens to be in a work tree, and the module source verbatim. A hash lets you verify a file you already have; the text is what answers “what produced this” when nobody can find the module a month later.

ROIs as segmentation layers¶

neu-mark roi-survey measures every ROI on a node before a spec is written, and again when the ROIs change:

neu-mark roi-survey --src dvid://dvid.example.org/abc123 \
    --compare dvid://dvid.example.org/def456 --spec rois.toml --like s3://my-bucket/seg \
    --out survey.json

It prints the containment tree recovered from the blocks alone: each ROI’s parent is the smallest larger ROI holding most of it, marked with its share inside when it sticks out (below --contained, 95%), and each parent with how much of it its children cover. Then the overlaps between ROIs that do not nest, flagged past --max-overlap as roi-layers would refuse them in one layer; empty and identical ROIs; ROIs reaching past --like’s extent. --compare reports each ROI whose blocks changed since another node, --spec checks a spec against this node, and --out writes the whole record, every overlapping pair included.

neu-mark roi-layers writes a node’s DVID ROIs as precomputed segmentations to view beside the segmentation, one volume per layer of a spec, each with segment properties naming its ROIs:

neu-mark roi-layers --src dvid://dvid.example.org/abc123 --spec rois.toml \
    --like s3://my-bucket/seg --dst s3://my-bucket/seg/rois

A voxel holds one label, so ROIs that nest (a medulla inside its optic lobe) go in different layers. The spec (TOML) declares:

layers: rois = [...] for a list of ROIs; tree = "top" and tree = "parts" for the top-level neuropils of its [neuropils] tree and their parts; from = LAYER for one label covering another layer. A parts layer is cut from the top layer’s regions, so the two nest exactly; with fill_gaps, blocks no part claims keep the top-level id. A one-label layer (from, or rois naming one ROI) can be cleaned up before it is written and meshed: open_blocks (an opening by a ball of that radius in ROI blocks, dropping specks and thin spurs), then close_blocks (a closing, bridging gaps and dents), then fill_holes = true (enclosed cavities filled).
[neuropils]: each top-level neuropil and the ROIs that make it up. One with parts is their union, plus its own ROI where the node has one.
[ids], shared by every layer, so saved states keep meaning the same segments after a rebuild; [names] for descriptions; [ignore] for ROIs left out on purpose. Every ROI on the node must be in a layer or ignored.

Within a layer, overlap settles ROIs drawn slightly over each other: split (the default) gives a contested block to the ROI it lies deepest inside, order to the one listed first, error refuses. ROIs overlapping by more than max_overlap (default 20%) of the smaller are nested rather than imprecise, and are refused.

An ROI is whole 32-voxel blocks of the node’s scale-0 grid; --like names the segmentation whose level-0 voxel size and origin that grid is. <dst>/rois.json records the node, the spec’s hash, the ids and each layer’s overlaps. Meshing is neu-morpho’s: neu-morpho run --stages index,mesh on each layer, with --occupancy-scale -1.

Layering¶

Sits above neu-vol (location for every write, dvid for node resolution, ops.provenance and ops.naming) and alongside neu-morpho, from which it imports nothing.