Match Seurat's violin: nrd0 bandwidth, trim, and points

[shanikawm]

Three defects in vln_plot. All checked against R 4.6.1 / Seurat 5.5.1
running locally, not against recollection of what Seurat does.

The bandwidth was 2.3x too wide

This is the real defect, and it is not the one I set out to fix.

The density used scipy's "scott" rule, which scales the sample standard
deviation. R's stats::density — and so geom_violin, and so Seurat — uses
bw.nrd0, which takes min(sd, IQR/1.34). Expression is zero-inflated, so the
IQR term is much the smaller of the two:

data	n	R nrd0	scipy scott	ratio
gamma	200	0.3901	0.4334	1.11x
gamma	2000	0.2680	0.3113	1.16x
zero-inflated	200	0.2069	0.4078	1.97x
zero-inflated	2000	0.1020	0.2503	2.45x

Over-smoothing by that much flattens the spike at zero, which is the shape of
single-cell expression. _bw_nrd0 now reproduces R to 5e-13 across seven
cases, including both branches of the zero-spread fallback chain.

bw.nrd0 divides the IQR by 1.34. The neighbouring rule bw.nrd, which R
also ships, uses 1.349 — and the first version of this fix used that one, a
silent 0.67% error wherever the IQR term wins. Diffing against R caught it.
That is why the constant is now commented in the source.

The violin was not trimmed

geom_violin(trim = TRUE) limits the density to the observed range. An
untrimmed gaussian KDE runs past it, so every violin tailed off below zero,
where expression cannot go. This was the originally-reported issue and it is the
smaller half.

Points were off by default

pt_size defaulted to 0. Seurat's VlnPlot passes pt.size = NULL, and
ExIPlot resolves it through AutoPointSize — min(1583/n, 1), read out of the
installed package. Points are now on by default and shrink as the cell count
grows (0.60 at pbmc3k's 2638 cells); pt_size=0 still omits them.

Validation

The drawn outline was compared against geom_violin(scale = "width", trim = TRUE)'s own computed polygon, pulled out of ggplot_build:

| group | max |Δ width| | support == data range |
|---|---|---|
| A | 1.4e-03 | yes |
| B | 6.4e-04 | yes |
| C | 2.1e-03 | yes |

0.2% of full width worst case; the residual is R's FFT-binned density() against
scipy's exact evaluation. Trim confirmed exactly.

One test earns its place specifically

Asserting _bw_nrd0 matches R says nothing about whether vln_plot calls
it. Substituting bw_method="scott" into the render path left all seventeen
other tests in the file green — the component was tested, the integration was
not.

The added test compares the drawn outline against both candidate densities and
requires the nrd0 one, with 590x separation between them (4e-4 against
2.3e-1). Full mutation sweep, all five caught:

Mutation	Result
scott instead of nrd0	caught (only by the new test)
divisor 1.34 → 1.349	caught, 4 tests
drop the trim	caught, 2 tests
area scaling instead of width	caught, 2 tests
points off by default	caught, 2 tests

Figures

Nine regenerated — exactly the nine drawn by the six generators that call
vln_plot, and no others. That the rest held still is also a check on #87: the
locked environment reproduces.

New: violin_width, and jitter_seed (default 0) so a committed figure
redraws identically.

Still deliberately different from Seurat: the median bar. geom_violin draws
none; this keeps drawing one.

1041 passed, 25 skipped.

🤖 Generated with Claude Code