Reproducing Pakhomov et al. (2023) from dandiset 000686¶
Paper. Pakhomov, Gudvangen, Mangalanathan, Kondratiev, Redondo & Semenov (2023), "Next generation
CANCAN focusing for remote stimulation by nanosecond electric pulses", Bioelectrochemistry 152,
108437, doi:10.1016/j.bioelechem.2023.108437.
The claims file (claims/000686/10.1016-j.bioelechem.2023.108437.yaml) was written from the accepted
author manuscript in PubMed Central (PMC10247516), Methods, Results, Discussion and legends; the paper
has no supplement and links no code.
Dataset. DANDI:000686, version 0.240515.1909, streamed from the archive. It is the only published version (the draft has the same 18 assets and was last modified a second after it was published), and the paper cites no version, so the most recent published version is pinned. The dandiset holds 18 NWB files, one per dish, each a single stitched YO-PRO-1 image of a cell monolayer inside the quadrupole array after treatment: the four groups of Fig. 5 (40 Packets and 40 Trains from e1, n = 4 each; 1P x 40 and 1T x 40 with rotation, n = 5 each), all Packet 2 or its train at 6.4 kV. No ROI measurements are deposited: the profiles are measured here from the images, with the ROI layout of Methods 2.4.
| Claim | Figure | Status here |
|---|---|---|
| fig3-field-simulation | 3 | blocked: data_missing (simulated field maps not deposited) |
| fig4-sham-no-uptake, fig4-unipolar-uptake, fig4-bipolar-weaker, fig4-8phase-cancellation | 4 | blocked: data_missing (single-cell time-lapse not deposited) |
| fig5a-trains-stronger | 5A | tested |
| fig5a-packet-center-enhanced | 5A | tested |
| fig5b-rotation-center-peak | 5B | tested |
| fig5c-normalized-contrast | 5C | tested |
| fig6a, fig6b, fig6-void-between, fig6-peak-narrower | 6 | blocked: data_missing (not in this dandiset; apparently in DANDI:000875) |
| fig7a-amplitude-scales-uniformly, fig7-amplitude-no-effect-on-distribution | 7 | blocked: data_missing (apparently in DANDI:001025) |
| fig8-center-equal-packets, fig8a-periphery-packets-1-3-higher, fig8b-normalized-center | 8 | blocked: data_missing (apparently in DANDI:001027) |
| fig9a-gap-reduces-periphery-more, fig9b-gap-boosts-targeting | 9 | blocked: data_missing (apparently in DANDI:001028) |
How the profiles are measured (helper cancan.py, which cites the source of every parameter).
Each image is streamed and reduced to 8 x 8 pixel block means. The four electrode footprints (the
"dark circles" of the file description) are found automatically. Along each diagonal, ROIs
run from the array centre (the mean of the four footprint centres) to each electrode: 56 per half
diagonal as in Methods 2.4 (224 ROIs along both diagonals), plus one on the centre, ending at
3.65/4.2 of the distance from the centre to the electrode's inner edge (Methods 2.4, +/- 3.65 mm; Results
3.5, 8.4-mm diagonal). Each ROI is one step (about 65 um) long along the diagonal. The paper gives
neither the ROI width across the diagonal nor how the electrode edge was placed, so six ROI geometries
are run: widths of 1, 3 and 5 steps, and the inner edge either at the footprint edge measured in the
image or at a fixed fraction of the centre-to-centre distance from the 1.5-mm needles (Methods 2.3).
The main tables use width 1 and the measured footprint; every claim's test is repeated under all six. The readings O1 to O6 of the claims file are run as it names them: the centre
as the centre ROI or a LOESS fit at 0 (spans 0.1, 0.2, 0.3), the periphery as the outermost ROI at
both ends or at the e1 end, raw values or values minus a background (the median inside the
footprints), normalization per dish or of the mean, and paired and unpaired t-tests where the
claims file asks for both. Tests are two-sided Student's t-tests (Methods 2.4) over dishes.
import collections
import re
from concurrent.futures import ThreadPoolExecutor
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import cancan
# ROI geometries (end of the ROIs, width across the diagonal in steps); MAIN is used in the main tables
GEOMS = [(e, w) for e in cancan.ROI_ENDS for w in cancan.ROI_WIDTHS]
MAIN = ("footprint", 1)
print("Dandiset", cancan.DANDISET_ID, "version", cancan.DANDISET_VERSION)
def load(r):
"""Stream one file, reduce it to 8 x 8 block means, find the electrodes and read the profiles."""
img, meta = cancan.load_image(r["url"])
el = cancan.find_electrodes(img)
h, w = (img.shape[0] // 8) * 8, (img.shape[1] // 8) * 8
coarse = np.nanmean(img[:h, :w].reshape(h // 8, 8, w // 8, 8), axis=(1, 3))
prof = {g: cancan.diagonal_profiles(img, el, across=g[1], end=g[0]) for g in GEOMS}
return dict(**r, meta=meta, el=el, prof={g: v["halves"] for g, v in prof.items()}, mm=prof[MAIN]["mm"],
bg=cancan.footprint_background(img, el), coarse=coarse)
import warnings
warnings.filterwarnings("ignore", "Mean of empty slice") # blocks outside the stitched tiles
with ThreadPoolExecutor(3) as ex:
D = list(ex.map(load, cancan.assets()))
print("loaded", len(D), "images")
Dandiset 000686 version 0.240515.1909
loaded 18 images
Data check¶
Counts by group and date, the metadata of the image arrays, the geometry found in each image, and a search for files deposited twice. Distances in mm assume 0.65 um per pixel (see the helper); no measure depends on it.
import collections
groups = ["40 Packets", "40 Trains", "1P x 40", "1T x 40"]
print("Dandiset", cancan.DANDISET_ID, "version", cancan.DANDISET_VERSION, "-", len(D), "NWB files")
print("Files per group:", dict(collections.Counter(d["group"] for d in D)))
print("Sessions by date:", dict(collections.Counter((d["group"], d["meta"]["session_start"][:10]) for d in D)))
m0 = D[0]["meta"]
print("\nImage array: shapes", {d["meta"]["shape"] for d in D}, "dtype", {d["meta"]["dtype"] for d in D},
"filters", {tuple(d["meta"]["filters"]) for d in D})
print("unit", {d["meta"]["unit"] for d in D}, "resolution", {d["meta"]["resolution"] for d in D},
"conversion", {d["meta"]["conversion"] for d in D}, "(no pixel size is given)")
print("species", {d["meta"]["species"] for d in D}, "(NCBITaxon_9913 is Bos taurus; the paper's cells are CHO-K1, hamster)")
print("Image description:", m0["image_description"])
print("\nTreatment named in each session description (pulses per packet or train, first-pulse kV, cycles):")
for g in groups:
descs = {re.sub(r"\s+", " ", d["meta"]["session_description"]) for d in D if d["group"] == g}
for s in descs:
print(f" {g}: {len(descs)} distinct;", re.findall(r"(?:packets|trains) of \d+ pulses|[\d.]+ kV|total of \d+ \w+", s))
rows = []
for d in D:
el, m = d["el"], d["meta"]
px_mm = cancan.BLOCK * cancan.UM_PER_PIXEL / 1000
rows.append(dict(group=d["group"], dish=d["dish"], start=m["session_start"][:16],
diag_centres_mm=el["diagonal"] * px_mm, footprint_diam_mm=2 * el["radius"] * px_mm,
diag_inner_edges_mm=(el["diagonal"] - 2 * el["radius"]) * px_mm,
square_err=el["square_error"], roi_end_mm=d["mm"][:, -1].mean(),
roi_end_needle_mm=np.mean([cancan.roi_positions(el, i, "needle")[1] for i in range(4)]) * cancan.N_HALF * px_mm,
roi_step_um=1000 * d["mm"][:, 1].mean(), background=d["bg"],
median=np.nanmedian(d["coarse"]), max=m["max"], frac_ge_65000=m["frac_ge_65000"],
frac_no_data=m["frac_zero"]))
T = pd.DataFrame(rows)
with pd.option_context("display.width", 250, "display.max_columns", 30, "display.precision", 3):
print(T.to_string(index=False))
print("\nRanges: centre-to-centre diagonal %.2f to %.2f mm, footprint diameter %.2f to %.2f mm, inner-edge diagonal "
"%.2f to %.2f mm (at %.2f um per pixel)" % (T.diag_centres_mm.min(), T.diag_centres_mm.max(), T.footprint_diam_mm.min(),
T.footprint_diam_mm.max(), T.diag_inner_edges_mm.min(), T.diag_inner_edges_mm.max(), cancan.UM_PER_PIXEL))
print("Tile geometry: 22323 px from 12 tiles is %.0f px per tile step, consistent with 2048-px ORCA-Flash4 tiles "
"at 10%% overlap (6.5 um / 10x = 0.65 um per pixel)." % ((22323 - 2048) / 11))
print("Profile values over all dishes and ROIs: %.0f to %.0f (raw block means, a.u.)" % (
min(np.nanmin(d["prof"][MAIN]) for d in D), max(np.nanmax(d["prof"][MAIN]) for d in D)))
print("NaN ROIs, all six ROI geometries:", sum(int(np.isnan(v).sum()) for d in D for v in d["prof"].values()))
print("Footprint radius against the fixed needle reading: ROI end with the measured footprint %.2f to %.2f mm, "
"with the needle reading %.2f to %.2f mm" % (T.roi_end_mm.min(), T.roi_end_mm.max(), T.roi_end_needle_mm.min(), T.roi_end_needle_mm.max()))
# Recordings deposited twice: correlation of the coarse images between every pair of files
C = np.array([d["coarse"].ravel() for d in D], dtype=float)
C = np.nan_to_num(C)
R = np.corrcoef(C)
iu = np.triu_indices(len(D), 1)
o = np.argsort(R[iu])[::-1][:3]
print("\nHighest image correlations between two files (1 would be a duplicate):")
for k in o:
i, j = iu[0][k], iu[1][k]
print(f" {D[i]['group']} {D[i]['dish']} / {D[j]['group']} {D[j]['dish']}: r = {R[i, j]:.3f}")
cent = np.array([d["el"]["centres"].ravel() for d in D])
dmin = min(np.abs(cent[i] - cent[j]).max() for i, j in zip(*iu))
print("Smallest largest-difference of footprint positions between two files: %.0f px" % dmin)
print("Subject ids are unique:", len({d['meta']['subject_id'] for d in D}) == len(D))
print("Subject age in the files:", {d["meta"]["age"] for d in D}, "(an age given for a cell line)")
from dandi.dandiapi import DandiAPIClient
with DandiAPIClient() as client:
summary = client.get_dandiset(cancan.DANDISET_ID, cancan.DANDISET_VERSION).get_raw_metadata()["assetsSummary"]
print("Dandiset assets summary: numberOfSubjects", summary.get("numberOfSubjects"), ", numberOfFiles", summary.get("numberOfFiles"),
"; distinct subject ids in the files:", len({d['meta']['subject_id'] for d in D}))
print("Timing: each file holds one image taken after treatment; there is no time series or trial table, "
"so there is no clock to check.")
Dandiset 000686 version 0.240515.1909 - 18 NWB files
Files per group: {'40 Packets': 4, '40 Trains': 4, '1P x 40': 5, '1T x 40': 5}
Sessions by date: {('40 Packets', '2021-08-03'): 2, ('40 Packets', '2021-08-06'): 2, ('40 Trains', '2021-08-03'): 2, ('40 Trains', '2021-08-06'): 2, ('1P x 40', '2021-08-06'): 3, ('1P x 40', '2021-08-03'): 2, ('1T x 40', '2021-08-06'): 3, ('1T x 40', '2021-08-03'): 2}
Image array: shapes {(1, 22323, 22323)} dtype {'uint16'} filters {('deflate',)}
unit {'n.a.'} resolution {-1.0} conversion {1.0} (no pixel size is given)
species {'http://purl.obolibrary.org/obo/NCBITaxon_9913'} (NCBITaxon_9913 is Bos taurus; the paper's cells are CHO-K1, hamster)
Image description: IX83 microscope used; MS-2000 scanning stage; X-Cite 110LED illuminator; ORCA-Flash4 sCMOS camera;CellSens software for automation; ex/em 480/535nm for detecting electroporation using 1uM YOPRO-1; dark circles indicate electrode footprints; 144 images captured in a 12x12 grid from a six-well plate; all images automatically stitched for a final composite; 10×, 0.38NA objective used for high-resolution imaging.
Treatment named in each session description (pulses per packet or train, first-pulse kV, cycles):
40 Packets: 1 distinct; ['packets of 15 pulses', '6.4 kV', 'total of 40 packets']
40 Trains: 1 distinct; ['trains of 4 pulses', '6.4 kV', 'total of 40 cycles']
1P x 40: 1 distinct; ['packets of 15 pulses', '6.4 kV', 'total of 10 cycles']
1T x 40: 1 distinct; ['trains of 4 pulses', '6.4 kV', 'total of 10 cycles']
group dish start diag_centres_mm footprint_diam_mm diag_inner_edges_mm square_err roi_end_mm roi_end_needle_mm roi_step_um background median max frac_ge_65000 frac_no_data
40 Packets dish09 2021-08-03T13:44 9.861 1.331 8.530 0.103 3.706 3.636 66.186 208.633 716.656 65443 4.469e-06 0.164
40 Packets dish10 2021-08-03T14:15 9.875 1.269 8.606 0.101 3.740 3.641 66.781 230.977 818.992 65443 7.467e-06 0.162
40 Packets dish10-2 2021-08-06T14:12 9.849 1.253 8.596 0.104 3.736 3.632 66.710 275.902 858.516 65448 2.583e-05 0.163
40 Packets dish11 2021-08-06T14:17 9.745 1.310 8.435 0.119 3.666 3.594 65.467 210.902 682.762 65443 7.371e-06 0.163
40 Trains dish11 2021-08-03T14:21 9.774 1.269 8.505 0.083 3.696 3.604 66.003 362.277 1312.753 65449 1.056e-03 0.162
40 Trains dish12 2021-08-03T14:27 9.869 1.378 8.491 0.088 3.690 3.639 65.897 261.699 2762.344 65449 6.608e-04 0.162
40 Trains dish14 2021-08-06T14:58 9.977 1.321 8.657 0.053 3.762 3.679 67.172 278.215 1417.550 65448 1.107e-03 0.162
40 Trains dish9 2021-08-06T13:39 9.865 1.404 8.461 0.050 3.677 3.637 65.659 295.266 1872.474 65449 5.512e-04 0.162
1P x 40 dish15 2021-08-06T15:04 9.820 1.248 8.572 0.013 3.725 3.621 66.517 211.836 746.929 65446 1.716e-05 0.163
1P x 40 dish2 2021-08-03T12:07 9.949 1.321 8.629 0.097 3.750 3.668 66.956 221.715 647.979 65442 7.469e-06 0.163
1P x 40 dish2-2 2021-08-06T12:04 9.960 1.149 8.810 0.061 3.829 3.672 68.367 207.547 638.642 65445 2.054e-05 0.163
1P x 40 dish4 2021-08-03T12:46 9.712 1.295 8.417 0.070 3.660 3.583 65.354 237.617 816.615 65444 1.196e-05 0.163
1P x 40 dish6 2021-08-06T12:56 9.866 1.269 8.597 0.068 3.736 3.638 66.712 181.531 618.406 65444 3.006e-06 0.163
1T x 40 dish19 2021-08-06T16:38 9.977 1.378 8.599 0.053 3.737 3.679 66.730 329.664 4077.574 65448 4.386e-04 0.161
1T x 40 dish3 2021-08-06T12:10 9.937 1.295 8.642 0.118 3.757 3.665 67.080 282.879 2677.009 65446 9.898e-05 0.161
1T x 40 dish5 2021-08-03T12:52 10.044 1.440 8.603 0.056 3.739 3.703 66.766 322.531 3689.872 65447 2.683e-04 0.161
1T x 40 dish5-2 2021-08-06T12:49 9.947 1.253 8.693 0.055 3.778 3.667 67.456 301.395 2550.873 65446 9.316e-05 0.161
1T x 40 dish7 2021-08-03T13:33 10.073 1.331 8.742 0.099 3.799 3.714 67.843 334.793 3856.193 65447 1.292e-04 0.161
Ranges: centre-to-centre diagonal 9.71 to 10.07 mm, footprint diameter 1.15 to 1.44 mm, inner-edge diagonal 8.42 to 8.81 mm (at 0.65 um per pixel)
Tile geometry: 22323 px from 12 tiles is 1843 px per tile step, consistent with 2048-px ORCA-Flash4 tiles at 10% overlap (6.5 um / 10x = 0.65 um per pixel).
Profile values over all dishes and ROIs: 403 to 24745 (raw block means, a.u.)
NaN ROIs, all six ROI geometries: 0
Footprint radius against the fixed needle reading: ROI end with the measured footprint 3.66 to 3.83 mm, with the needle reading 3.58 to 3.71 mm
Highest image correlations between two files (1 would be a duplicate):
40 Trains dish14 / 40 Trains dish9: r = 0.815
40 Trains dish12 / 40 Trains dish9: r = 0.745
40 Trains dish11 / 40 Trains dish14: r = 0.734
Smallest largest-difference of footprint positions between two files: 35 px
Subject ids are unique: True
Subject age in the files: {'P1D'} (an age given for a cell line)
Dandiset assets summary: numberOfSubjects 36 , numberOfFiles 18 ; distinct subject ids in the files: 18 Timing: each file holds one image taken after treatment; there is no time series or trial table, so there is no clock to check.
Images¶
Each dish's image at 1/64 of the resolution, with the four footprints found and the ROIs of both diagonals (every 4th ROI drawn). Brightness is scaled to each image's 99.5th percentile.
fig, axs = plt.subplots(3, 6, figsize=(24, 12.5))
for a, d in zip(axs.flat, D):
x = d["coarse"]
a.imshow(x, cmap="gray", vmax=np.nanpercentile(x, 99.5))
for i in range(4):
pts, _ = cancan.roi_positions(d["el"], i)
a.plot(pts[::4, 1] / 8, pts[::4, 0] / 8, ".", ms=1.5, color="C1" if i in (1, 3) else "C2")
c = d["el"]["centres"] / 8
a.plot(c[:, 1], c[:, 0], "o", mfc="none", color="C3", ms=12)
a.set_title(f"{d['group']} {d['dish']}", fontsize=10); a.axis("off")
plt.tight_layout(); plt.show()
Which electrode is e1¶
The files do not say which footprint is e1. Fig. 5A and Results 3.2 need it (the profile runs from e1 to e3). Two readings (i) the electrode with the brightest outermost ROIs in the non-rotated dishes; (ii) the top-right electrode of the image, the mirror image of the top-left electrode energized in Fig. 3 (an inverted microscope sees the dish from below). Reading (i) uses the uptake near the electrode and so is partly circular for the train profile, which is why reading (ii) is shown with it.
# Which electrode is e1. Reading (i): e1 is the electrode with the brightest outer ROIs in the
# non-rotated dishes (for trains, Results 3.2 says uptake follows the field, strongest near e1).
# Reading (ii): e1 is the top-right electrode of the image, the mirror of the top-left electrode
# energized in Fig. 3, as an inverted microscope sees the dish from below.
names = ["top left", "top right", "bottom right", "bottom left"]
print("Mean of the 3 outermost ROIs at each electrode (image clockwise from top left):")
bright = []
for d in D:
if d["rotated"]:
continue
v = d["prof"][MAIN][:, -3:].mean(1)
bright.append(int(np.argmax(v)))
print(f" {d['group']:10s} {d['dish']:9s}", np.round(v).astype(int), "brightest:", names[bright[-1]])
print("Brightest electrode over the 8 non-rotated dishes:", dict(collections.Counter(names[b] for b in bright)))
E1, E3 = 1, 3 # top right and bottom left; both readings agree when every dish is top right
E2, E4 = 2, 0 # the other diagonal; E2 next to E1 clockwise in the image (a labelling choice)
print("e1 = top right under both readings:", all(b == E1 for b in bright))
Mean of the 3 outermost ROIs at each electrode (image clockwise from top left):
40 Packets dish09 [1185 3429 1055 803] brightest: top right
40 Packets dish10 [2233 8677 2360 1354] brightest: top right
40 Packets dish10-2 [2029 5397 1684 1735] brightest: top right
40 Packets dish11 [1003 2855 1136 1173] brightest: top right
40 Trains dish11 [ 9428 16987 11850 1101] brightest: top right
40 Trains dish12 [ 6610 13804 8553 1190] brightest: top right
40 Trains dish14 [ 7204 13156 11509 1250] brightest: top right
40 Trains dish9 [ 5203 13467 7966 1130] brightest: top right
Brightest electrode over the 8 non-rotated dishes: {'top right': 8}
e1 = top right under both readings: True
Mean profiles¶
Fig. 5A (one diagonal, e3 at left, e1 at right) and Fig. 5B (both diagonals averaged by signed distance). Mean +/- s.e. over dishes.
O1 = ["roi"] + [f"loess{s}" for s in cancan.LOESS_SPANS]
O4 = {"raw": False, "bg": True}
def group(g):
return [d for d in D if d["group"] == g]
def prof(d, bg=False, geom=MAIN):
"""Signed profile along e1-e3 (index 0 at e3, 56 at centre, 112 at e1)."""
return cancan.signed_profile(d["prof"][geom], E1, E3) - (d["bg"] if bg else 0.0)
def prof_both(d, bg=False, geom=MAIN):
"""Profile averaged over both diagonals by signed distance (Fig. 5B legend); the positive end
averages the e1 and e2 ends, the negative end the e3 and e4 ends."""
h = d["prof"][geom]
return (cancan.signed_profile(h, E1, E3) + cancan.signed_profile(h, E2, E4)) / 2 - (d["bg"] if bg else 0.0)
def trough_index(S, between, loo):
"""The trough ROI: the lowest point of the group mean profile in `between`. With `loo`, each
dish's trough is found on the mean of the other dishes (leave one out), so that a dish's own
noise does not pick where it is measured."""
if not loo:
return np.full(len(S), between[np.argmin(S.mean(0)[between])])
return np.array([between[np.argmin(np.delete(S, i, 0).mean(0)[between])] for i in range(len(S))])
def test_5a_trains(geom=MAIN):
P, Tr = group("40 Packets"), group("40 Trains")
rows = []
for o1 in O1:
for o4, bg in O4.items():
vals = {}
for name, G in [("P", P), ("T", Tr)]:
Y = [prof(d, bg, geom) for d in G]
vals[name] = dict(e1=np.array([y[-1] for y in Y]), centre=np.array([cancan.centre_value(y, o1) for y in Y]),
e3=np.array([y[0] for y in Y]))
r = dict(O1=o1, O4=o4)
ok = True
for place in ["e1", "centre", "e3"]:
t, p = cancan.ttest(vals["T"][place], vals["P"][place])
r[f"T/P {place}"] = f"{vals['T'][place].mean():.0f}/{vals['P'][place].mean():.0f} p={p:.3g}"
r[f"ok {place}"] = (t > 0) and (p < 0.05)
r[f"p {place}"] = p
ok &= r[f"ok {place}"]
t, p = cancan.ttest(vals["T"]["e1"], vals["T"]["centre"])
tp, pp = cancan.ttest(vals["T"]["e1"], vals["T"]["centre"], paired=True)
r["T e1>centre unpaired"] = f"t={t:.2f} p={p:.3g}"
r["(paired variant)"] = f"t={tp:.2f} p={pp:.3g}"
r["ok fall"] = (t > 0) and (p < 0.05)
r["e3 bar gap (s.e.)"] = round(cancan.errorbar_gap(vals["T"]["e3"], vals["P"]["e3"])[0], 2)
r["holds"] = ok and r["ok fall"]
rows.append(r)
return pd.DataFrame(rows)
def test_5a_packet(geom=MAIN, loo=False):
P = group("40 Packets")
rows = []
for o1 in O1:
for o4, bg in O4.items():
Y = np.array([prof(d, bg, geom) for d in P])
S = np.array([cancan.smooth(y, o1) for y in Y])
k = trough_index(S, np.arange(57, 112), loo) # strictly between the centre (56) and the e1 end (112)
centre, trough, e1 = S[:, 56], S[np.arange(len(S)), k], Y[:, -1]
for paired in (False, True):
t1, p1 = cancan.ttest(centre, trough, paired)
t2, p2 = cancan.ttest(e1, centre, paired)
rows.append(dict(O1=o1, O4=o4, test="paired" if paired else "unpaired",
trough_mm=round(float(x_mm[k[0]]), 2), centre=round(centre.mean()), trough=round(trough.mean()),
e1_end=round(e1.mean()), **{"centre>trough": f"t={t1:.2f} p={p1:.3g}",
"e1>centre": f"t={t2:.2f} p={p2:.3g}"}, p_centre=p1, p_e1=p2,
holds_centre=(t1 > 0) and (p1 < 0.05), holds_e1=(t2 > 0) and (p2 < 0.05),
gap_centre_trough=round(cancan.errorbar_gap(centre, trough)[0], 2),
gap_e1_centre=round(cancan.errorbar_gap(e1, centre)[0], 2)))
R = pd.DataFrame(rows)
R["holds"] = R.holds_centre & R.holds_e1
return R
def test_5b(geom=MAIN, loo=False):
RP, RT = group("1P x 40"), group("1T x 40")
rows = []
for o1 in O1:
for o4, bg in O4.items():
YP = np.array([prof_both(d, bg, geom) for d in RP]); SP = np.array([cancan.smooth(y, o1) for y in YP])
YT = np.array([prof_both(d, bg, geom) for d in RT]); ST = np.array([cancan.smooth(y, o1) for y in YT])
k = trough_index(SP, np.r_[1:56, 57:112], loo) # strictly between the centre and either end
cP, trough = SP[:, 56], SP[np.arange(len(SP)), k]
cT = ST[:, 56]
for o2 in ["both", "e1+e2"]:
perT = (YT[:, 0] + YT[:, -1]) / 2 if o2 == "both" else YT[:, -1]
for paired in (False, True):
t1, p1 = cancan.ttest(cP, trough, paired)
t2, p2 = cancan.ttest(perT, cT, paired)
rows.append(dict(O1=o1, O4=o4, O2=o2, test="paired" if paired else "unpaired",
P_centre=round(cP.mean()), P_trough=round(trough.mean()), trough_mm=round(float(x_mm[k[0]]), 2),
**{"P centre>trough": f"t={t1:.2f} p={p1:.3g}"},
T_centre=round(cT.mean()), T_periph=round(perT.mean()),
**{"T periph>centre": f"t={t2:.2f} p={p2:.3g}"}, p_P=p1, p_T=p2,
holds_P=(t1 > 0) and (p1 < 0.05), holds_T=(t2 > 0) and (p2 < 0.05),
T_centre_lowest=bool(np.argmin(ST.mean(0)) == 56),
gap_P=round(cancan.errorbar_gap(cP, trough)[0], 2),
gap_T=round(cancan.errorbar_gap(perT, cT)[0], 2)))
R = pd.DataFrame(rows)
R["holds"] = R.holds_P & R.holds_T
return R
def test_5c(geom=MAIN, keep=None):
RP, RT = group("1P x 40"), group("1T x 40")
rows = []
for o1 in O1:
for o4, bg in O4.items():
for o2 in ["both", "e1+e2"]:
per, cen = {}, {}
for name, G in [("P", RP), ("T", RT)]:
Y = np.array([prof_both(d, bg, geom) for d in G])
cen[name] = np.array([cancan.centre_value(y, o1) for y in Y])
per[name] = (Y[:, 0] + Y[:, -1]) / 2 if o2 == "both" else Y[:, -1]
if keep is not None and o1 == "roi" and o2 == "both":
keep[(o4, name)] = 100 * Y / per[name][:, None]
for o3 in ["a", "b"]:
if o3 == "a":
nP, nT = 100 * cen["P"] / per["P"], 100 * cen["T"] / per["T"]
vP, vT = nP.mean(), nT.mean()
else:
nP, nT = 100 * cen["P"] / per["P"].mean(), 100 * cen["T"] / per["T"].mean()
vP, vT = 100 * cen["P"].mean() / per["P"].mean(), 100 * cen["T"].mean() / per["T"].mean()
t, p = cancan.ttest(nP, nT)
rows.append(dict(O1=o1, O4=o4, O2=o2, O3=o3, P_norm_centre=round(vP, 1), T_norm_centre=round(vT, 1),
diff_points=round(vP - vT, 1), t=round(t, 2), p=float(f"{p:.3g}"),
holds=(t > 0) and (p < 0.05), P_above_100=vP > 100,
gap=round(cancan.errorbar_gap(nP, nT)[0], 2)))
return pd.DataFrame(rows)
def by_geometry(test, cols, **kw):
"""Run one claim's test under every ROI geometry and count, for each, the readings where each
condition holds (of all readings of the claims file)."""
rows = []
for g in GEOMS:
R = test(g, **kw)
rows.append(dict(end=g[0], width=g[1], readings=len(R), **{c: int(R[c].sum()) for c in cols}))
return pd.DataFrame(rows)
x_mm = np.concatenate([-D[0]["mm"][E3][::-1], D[0]["mm"][E1][1:]]) # distances of the first dish, for labels
fig, ax = plt.subplots(1, 2, figsize=(13, 4.5))
for g, c in [("40 Packets", "C0"), ("40 Trains", "C3")]:
Y = np.array([prof(d) for d in group(g)])
xs = np.mean([np.concatenate([-d["mm"][E3][::-1], d["mm"][E1][1:]]) for d in group(g)], 0)
ax[0].errorbar(xs, Y.mean(0), Y.std(0, ddof=1) / np.sqrt(len(Y)), color=c, lw=1, label=f"{g} (n={len(Y)})")
ax[0].plot(xs, cancan.loess(Y.mean(0), 0.2), color="k", lw=0.8)
print(f"{g}: mean profile range {Y.mean(0).min():.0f} to {Y.mean(0).max():.0f} a.u.")
ax[0].set(xlabel="distance from centre toward e1 (mm, at 0.65 um/px)", ylabel="YP intensity (a.u.)", title="Fig. 5A: e3 (left) to e1 (right)")
ax[0].legend()
for g, c in [("1P x 40", "C0"), ("1T x 40", "C3")]:
Y = np.array([prof_both(d) for d in group(g)])
ax[1].errorbar(x_mm, Y.mean(0), Y.std(0, ddof=1) / np.sqrt(len(Y)), color=c, lw=1, label=f"{g} (n={len(Y)})")
ax[1].plot(x_mm, cancan.loess(Y.mean(0), 0.2), color="k", lw=0.8)
print(f"{g}: mean profile range {Y.mean(0).min():.0f} to {Y.mean(0).max():.0f} a.u.")
ax[1].set(xlabel="distance from centre (mm)", title="Fig. 5B: averaged over both diagonals")
ax[1].legend(); plt.tight_layout(); plt.show()
print("All points are inside the axes (autoscaled). Black lines: LOESS, span 0.2, of the mean profile.")
40 Packets: mean profile range 690 to 6339 a.u. 40 Trains: mean profile range 789 to 19389 a.u. 1P x 40: mean profile range 1032 to 3581 a.u. 1T x 40: mean profile range 4884 to 11137 a.u.
All points are inside the axes (autoscaled). Black lines: LOESS, span 0.2, of the mean profile.
Claims tested¶
fig5a-trains-stronger¶
Without rotation, the control of 40 unipolar trains from e1 causes much stronger electroporation throughout the array than 40 Packet 2 packets started from e1. The uptake after trains follows the field strength distribution and falls with distance from e1.
Criterion (claims file): trains above packets at the e1 end (outermost ROI), the centre (O1) and the e3 end, each by an unpaired two-sided Student's t-test at p < 0.05, and the train e1-end value above the train centre by the same test; under O1 (a, b with spans 0.1, 0.2, 0.3) and O4 (raw, minus footprint background).
The test is run under every reading in the main ROI geometry (width 1, measured footprint) and then under all six ROI geometries.
P, Tr = group("40 Packets"), group("40 Trains")
R5a = test_5a_trains()
with pd.option_context("display.width", 250, "display.max_columns", 20):
print(R5a.drop(columns=["ok e1", "ok centre", "ok e3", "ok fall", "p e1", "p centre", "p e3"]).to_string(index=False))
print("\nReadings where every condition holds:", int(R5a.holds.sum()), "of", len(R5a))
print("\nAll readings under every ROI geometry: readings where each condition holds")
print(by_geometry(test_5a_trains, ["ok e1", "ok centre", "ok e3", "ok fall", "holds"]).to_string(index=False))
for g in GEOMS:
R = test_5a_trains(g)
print(f" {g}: e3 end T/P", R["T/P e3"].iloc[0], "(raw);", R["T/P e3"].iloc[1], "(bg); e1 end", R["T/P e1"].iloc[0],
"(raw); p over all 8 readings: e1 end %.3g to %.3g, centre %.3g to %.3g, e3 end %.3g to %.3g" % (
R["p e1"].min(), R["p e1"].max(), R["p centre"].min(), R["p centre"].max(), R["p e3"].min(), R["p e3"].max()))
print("\nPer dish, raw ROI values (e1 end, centre, e3 end) at widths 1, 3 and 5 (measured footprint):")
for name, G in [("40 Packets", P), ("40 Trains", Tr)]:
for d in G:
print(f" {name:10s} {d['dish']:9s}", " ".join("%6.0f %6.0f %6.0f" % tuple(prof(d, geom=("footprint", w))[[-1, 56, 0]])
for w in cancan.ROI_WIDTHS))
# Train profile shape: ROI means of the group mean profile from the e3 end to the e1 end, every 8 ROIs
Ym = np.mean([prof(d) for d in Tr], 0)
print("Train mean profile every 8th ROI from e3 (-) to e1 (+):", np.round(Ym[::8]).astype(int))
print("Trains against packets are compared between different dishes, so only the unpaired test applies (claims file, O6).")
O1 O4 T/P e1 T/P centre T/P e3 T e1>centre unpaired (paired variant) e3 bar gap (s.e.) holds
roi raw 13282/6339 p=0.0831 5998/2576 p=0.0657 1138/1299 p=0.558 t=2.57 p=0.0426 t=4.27 p=0.0236 -0.76 False
roi bg 12982/6108 p=0.0843 5699/2345 p=0.069 839/1068 p=0.373 t=2.59 p=0.0413 t=4.27 p=0.0236 -0.35 False
loess0.1 raw 13282/6339 p=0.0831 6385/2905 p=0.0146 1138/1299 p=0.558 t=2.56 p=0.0428 t=3.15 p=0.0514 -0.76 False
loess0.1 bg 12982/6108 p=0.0843 6086/2673 p=0.0156 839/1068 p=0.373 t=2.58 p=0.0416 t=3.15 p=0.0514 -0.35 False
loess0.2 raw 13282/6339 p=0.0831 6193/2731 p=0.00534 1138/1299 p=0.558 t=2.70 p=0.0355 t=3.25 p=0.0474 -0.76 False
loess0.2 bg 12982/6108 p=0.0843 5894/2499 p=0.00567 839/1068 p=0.373 t=2.72 p=0.0344 t=3.25 p=0.0474 -0.35 False
loess0.3 raw 13282/6339 p=0.0831 6293/2520 p=0.00332 1138/1299 p=0.558 t=2.67 p=0.0372 t=3.22 p=0.0485 -0.76 False
loess0.3 bg 12982/6108 p=0.0843 5993/2288 p=0.00348 839/1068 p=0.373 t=2.69 p=0.0361 t=3.22 p=0.0485 -0.35 False
Readings where every condition holds: 0 of 8
All readings under every ROI geometry: readings where each condition holds
end width readings ok e1 ok centre ok e3 ok fall holds
footprint 1 8 0 6 0 8 0
footprint 3 8 8 8 0 8 0
footprint 5 8 8 8 0 8 0
needle 1 8 8 6 0 8 0
needle 3 8 8 8 0 8 0
needle 5 8 8 8 0 8 0
('footprint', 1): e3 end T/P 1138/1299 p=0.558 (raw); 839/1068 p=0.373 (bg); e1 end 13282/6339 p=0.0831 (raw); p over all 8 readings: e1 end 0.0831 to 0.0843, centre 0.00332 to 0.069, e3 end 0.373 to 0.558
('footprint', 3): e3 end T/P 1189/1385 p=0.522 (raw); 889/1153 p=0.367 (bg); e1 end 12913/5947 p=0.0244 (raw); p over all 8 readings: e1 end 0.0244 to 0.0247, centre 0.00402 to 0.0237, e3 end 0.367 to 0.522
('footprint', 5): e3 end T/P 1209/1327 p=0.644 (raw); 910/1096 p=0.437 (bg); e1 end 13087/5903 p=0.0153 (raw); p over all 8 readings: e1 end 0.0153 to 0.0154, centre 0.00607 to 0.0215, e3 end 0.437 to 0.644
('needle', 1): e3 end T/P 1291/1180 p=0.638 (raw); 992/948 p=0.858 (bg); e1 end 13576/4763 p=0.000702 (raw); p over all 8 readings: e1 end 0.000702 to 0.000719, centre 0.00321 to 0.066, e3 end 0.638 to 0.858
('needle', 3): e3 end T/P 1176/1287 p=0.508 (raw); 877/1056 p=0.268 (bg); e1 end 14889/4757 p=0.000252 (raw); p over all 8 readings: e1 end 0.000246 to 0.000252, centre 0.00394 to 0.0211, e3 end 0.268 to 0.508
('needle', 5): e3 end T/P 1236/1301 p=0.753 (raw); 936/1069 p=0.496 (bg); e1 end 14849/4870 p=0.000405 (raw); p over all 8 readings: e1 end 0.000389 to 0.000405, centre 0.00602 to 0.0254, e3 end 0.496 to 0.753
Per dish, raw ROI values (e1 end, centre, e3 end) at widths 1, 3 and 5 (measured footprint):
40 Packets dish09 4593 3941 897 4879 3270 892 4524 3220 971
40 Packets dish10 12931 3648 1282 10856 3400 1270 10758 3409 1224
40 Packets dish10-2 4601 1091 1971 4993 2124 2085 5149 2716 1906
40 Packets dish11 3233 1627 1048 3058 2141 1292 3181 2467 1210
40 Trains dish11 20757 8903 1370 17468 8547 1345 17090 7701 1391
40 Trains dish12 11154 7506 993 12821 9381 764 12003 9130 802
40 Trains dish14 10245 2994 938 10651 4554 1371 11576 4142 1352
40 Trains dish9 10972 4590 1252 10713 4617 1275 11679 5305 1292
Train mean profile every 8th ROI from e3 (-) to e1 (+): [ 1138 789 1029 1913 2298 2136 5821 5998 8447 12400 15721 17120
14168 17953 13282]
Trains against packets are compared between different dishes, so only the unpaired test applies (claims file, O6).
Result: not reproduced. The criterion fails under all 8 readings in every one of the six ROI geometries, because at the e3 end the trains do not exceed the packets in any of them (main geometry 1138 against 1299 a.u., p = 0.558; 839 against 1068 after background, p = 0.373; other geometries p = 0.268 to 0.858). In the main geometry the trains exceed the packets in the centre under the 6 LOESS readings (p = 0.00332 to 0.0156), not with the centre ROI (p = 0.0657 and 0.069), and near e1 only in the mean (13282 against 6339 a.u., p = 0.0831); with wider ROIs or the needle reading of the ROI end, the e1-end difference is significant in all 8 readings (p = 0.000246 to 0.0247 over the five other geometries) and the centre one in 6 or 8. The train e1-end value exceeds the train centre in all 8 main readings (unpaired p = 0.0344 to 0.0428) and in every geometry, and the train mean profile falls from 13282 to 17953 a.u. near e1 to 789 to 1138 a.u. near e3. So the second sentence of the statement holds and the first holds in the centre and near e1, not at the far end of the diagonal.
fig5a-packet-center-enhanced¶
Without rotation, Packet 2 enhances electroporation in the center of the array despite the weakening of the field with distance from e1. Uptake near e1 still stays stronger than in the center.
Criterion (claims file): on the 40-packet e1-e3 profile, the centre (O1) above the lowest value between the centre and the e1 end, and the e1-end value above the centre, each by a two-sided t-test at p < 0.05 over dishes, paired and unpaired (O6), under all readings of O1 and O4. The trough is the lowest point of the group's mean profile (ROI values or the mean of the dishes' LOESS fits, as O1) strictly between the centre and the e1 end; each dish's value at that ROI enters the test.
Variant (verifier, round 1): the trough chosen on the same dishes it is tested on favours centre > trough. A leave-one-out variant finds each dish's trough on the mean of the other dishes.
R5ap = test_5a_packet()
with pd.option_context("display.width", 250, "display.max_columns", 20):
print(R5ap.drop(columns=["p_centre", "p_e1"]).to_string(index=False))
print("\nReadings where centre > trough holds:", int(R5ap.holds_centre.sum()), "of", len(R5ap),
"; e1 end > centre holds:", int(R5ap.holds_e1.sum()), "; both:", int(R5ap.holds.sum()))
print("\nAll readings under every ROI geometry (trough on the group mean, and leave one out):")
G5ap = by_geometry(test_5a_packet, ["holds_centre", "holds_e1", "holds"])
L5ap = by_geometry(lambda g: test_5a_packet(g, loo=True), ["holds_centre"]).rename(columns={"holds_centre": "holds_centre_loo"})
G5ap["holds_centre_loo"] = L5ap.holds_centre_loo
G5ap["p e1>centre"] = [f"{test_5a_packet(g).p_e1.min():.3g} to {test_5a_packet(g).p_e1.max():.3g}" for g in GEOMS]
print(G5ap.to_string(index=False))
print("Per dish (raw ROI): centre, e1 end, ratio e1/centre")
for d in P:
y = prof(d)
print(f" {d['dish']:9s} {y[56]:7.0f} {y[-1]:7.0f} {y[-1] / y[56]:5.2f}")
Y = np.array([prof(d) for d in P])
print("Packet mean profile every 8th ROI from e3 (-) to e1 (+):", np.round(Y.mean(0)[::8]).astype(int))
O1 O4 test trough_mm centre trough e1_end centre>trough e1>centre holds_centre holds_e1 gap_centre_trough gap_e1_centre holds
roi raw unpaired 1.67 2576 1626 6339 t=1.26 p=0.255 t=1.61 p=0.158 False False -0.01 0.37 False
roi raw paired 1.67 2576 1626 6339 t=1.17 p=0.325 t=1.95 p=0.147 False False -0.01 0.37 False
roi bg unpaired 1.67 2345 1395 6108 t=1.24 p=0.261 t=1.61 p=0.158 False False -0.03 0.37 False
roi bg paired 1.67 2345 1395 6108 t=1.17 p=0.325 t=1.95 p=0.147 False False -0.03 0.37 False
loess0.1 raw unpaired 1.60 2905 1831 6339 t=4.46 p=0.00428 t=1.54 p=0.175 True False 3.37 0.44 False
loess0.1 raw paired 1.60 2905 1831 6339 t=5.14 p=0.0143 t=1.66 p=0.195 True False 3.37 0.44 False
loess0.1 bg unpaired 1.60 2673 1600 6108 t=4.37 p=0.00474 t=1.54 p=0.175 True False 3.29 0.44 False
loess0.1 bg paired 1.60 2673 1600 6108 t=5.14 p=0.0143 t=1.66 p=0.195 True False 3.29 0.44 False
loess0.2 raw unpaired 1.80 2731 1911 6339 t=4.33 p=0.00493 t=1.62 p=0.156 True False 3.48 0.55 False
loess0.2 raw paired 1.80 2731 1911 6339 t=3.64 p=0.0357 t=1.72 p=0.184 True False 3.48 0.55 False
loess0.2 bg unpaired 1.80 2499 1680 6108 t=4.35 p=0.0048 t=1.62 p=0.156 True False 3.53 0.55 False
loess0.2 bg paired 1.80 2499 1680 6108 t=3.64 p=0.0357 t=1.72 p=0.184 True False 3.53 0.55 False
loess0.3 raw unpaired 1.74 2520 1978 6339 t=3.46 p=0.0135 t=1.72 p=0.137 True False 2.48 0.66 False
loess0.3 raw paired 1.74 2520 1978 6339 t=2.96 p=0.0598 t=1.79 p=0.171 False False 2.48 0.66 False
loess0.3 bg unpaired 1.74 2288 1747 6108 t=3.51 p=0.0127 t=1.72 p=0.137 True False 2.60 0.66 False
loess0.3 bg paired 1.74 2288 1747 6108 t=2.96 p=0.0598 t=1.79 p=0.171 False False 2.60 0.66 False
Readings where centre > trough holds: 10 of 16 ; e1 end > centre holds: 0 ; both: 0
All readings under every ROI geometry (trough on the group mean, and leave one out):
end width readings holds_centre holds_e1 holds holds_centre_loo p e1>centre footprint 1 16 10 0 0 10 0.137 to 0.195 footprint 3 16 12 0 0 10 0.0955 to 0.155 footprint 5 16 14 0 0 14 0.098 to 0.164 needle 1 16 10 0 0 10 0.113 to 0.226 needle 3 16 10 0 0 10 0.135 to 0.231 needle 5 16 14 0 0 12 0.0533 to 0.116 Per dish (raw ROI): centre, e1 end, ratio e1/centre dish09 3941 4593 1.17 dish10 3648 12931 3.54 dish10-2 1091 4601 4.22 dish11 1627 3233 1.99 Packet mean profile every 8th ROI from e3 (-) to e1 (+): [1299 1764 1051 937 718 1877 1543 2576 2208 2241 1786 2094 2359 3559 6339]
Result: not reproduced. The e1-end value of the 40-packet profile exceeds the centre in the mean (6339 against 2576 a.u., raw ROIs) and in all 4 dishes (ratios 1.17 to 4.22), but under none of the 16 readings is the difference significant (p = 0.137 to 0.195 in the main geometry, 0.0533 to 0.231 over all six), and by the paper's error-bar rule it is not either (gap 0.37 to 0.66 s.e.). The centre exceeds the trough between centre and e1 under 10 of 16 main readings (every LOESS reading unpaired, spans 0.1 and 0.2 paired; p = 0.00428 to 0.0357), not with the centre ROI (p = 0.255 to 0.325) or span 0.3 paired (p = 0.0598); over the six geometries in 10 to 14 of 16, and in 10 to 14 with the trough found leave one out. As the criterion needs both comparisons and the second fails under every reading and geometry, the outcome is the same under all of them.
fig5b-rotation-center-peak¶
With rotation (1P x 40), Packet 2 smears the effect at the perimeter and brings out a remote effect in the center. The same rotation of unipolar trains (1T x 40) leaves a void of electroporation in the center.
Criterion (claims file): profiles averaged across both diagonals; the 1P x 40 centre (O1) above its lowest value between centre and periphery (on either side, found on the mean profile as above), and the 1T x 40 centre below its most peripheral value (O2: both ends, or the e1 end), each by a two-sided t-test at p < 0.05 (O6: paired and unpaired both run), under all readings of O1, O2 and O4. O2 (b) is the outermost ROI of the positive end of the averaged profile, which averages the e1 and e2 ends (a labelling choice: with rotation no end is e1's alone). The trough is also found leave one out, as for fig5a-packet-center-enhanced.
RP, RT = group("1P x 40"), group("1T x 40")
R5b = test_5b()
with pd.option_context("display.width", 260, "display.max_columns", 30):
print(R5b.drop(columns=["p_P", "p_T"]).to_string(index=False))
print("\n1P x 40 centre > trough holds in", int(R5b.holds_P.sum()), "of", len(R5b), "readings; 1T x 40 periphery > centre in",
int(R5b.holds_T.sum()), "; both in", int(R5b.holds.sum()))
print("1T x 40 centre is the lowest point of its mean profile (reported, not required):",
dict(collections.Counter(zip(R5b.O1, R5b.T_centre_lowest))))
for o4 in O4:
r = R5b[R5b.O4 == o4]
print(f"1T x 40, O4 {o4}: centre {r.T_centre.min()} to {r.T_centre.max()} a.u., most peripheral value {r.T_periph.min()} to {r.T_periph.max()} a.u. (main geometry)")
YT = np.array([prof_both(d) for d in RT])
print("1T x 40: mean profile minimum at %.2f mm (%.0f a.u.), centre %.0f" % (x_mm[np.argmin(YT.mean(0))], YT.mean(0).min(), YT.mean(0)[56]))
print("\nAll readings under every ROI geometry, split by paired and unpaired test (O6):")
rows = []
for g in GEOMS:
for loo in (False, True):
R = test_5b(g, loo=loo)
for test, r in R.groupby("test"):
rows.append(dict(end=g[0], width=g[1], trough="leave one out" if loo else "group mean", test=test, readings=len(r),
holds_P=int(r.holds_P.sum()), holds_T=int(r.holds_T.sum()), holds=int(r.holds.sum()),
max_p_P=float(f"{r.p_P.max():.3g}"), max_p_T=float(f"{r.p_T.max():.3g}"),
directions_as_claimed=bool(((r.P_centre > r.P_trough) & (r.T_periph > r.T_centre)).all())))
G5b = pd.DataFrame(rows)
with pd.option_context("display.width", 250, "display.max_columns", 20):
print(G5b.to_string(index=False))
print("Geometries where both comparisons hold under every reading (paired and unpaired):",
[g for g in GEOMS if G5b[(G5b.end == g[0]) & (G5b.width == g[1]) & (G5b.trough == "group mean")].holds.sum() == 32])
print("\nPer dish (raw ROI, averaged diagonals): centre, mean of both ends, at widths 1, 3 and 5 (measured footprint)")
for name, G in [("1P x 40", RP), ("1T x 40", RT)]:
for d in G:
print(f" {name:8s} {d['dish']:9s}", " ".join("%6.0f %6.0f" % (y[56], (y[0] + y[-1]) / 2)
for y in [prof_both(d, geom=("footprint", w)) for w in cancan.ROI_WIDTHS]))
O1 O4 O2 test P_centre P_trough trough_mm P centre>trough T_centre T_periph T periph>centre holds_P holds_T T_centre_lowest gap_P gap_T holds
roi raw both unpaired 3581 1032 -2.58 t=2.14 p=0.0647 4884 10910 t=2.98 p=0.0176 False True True 1.05 1.85 False
roi raw both paired 3581 1032 -2.58 t=2.02 p=0.114 4884 10910 t=2.81 p=0.0482 False True True 1.05 1.85 False
roi raw e1+e2 unpaired 3581 1032 -2.58 t=2.14 p=0.0647 4884 10989 t=3.10 p=0.0147 False True True 1.05 1.96 False
roi raw e1+e2 paired 3581 1032 -2.58 t=2.02 p=0.114 4884 10989 t=2.87 p=0.0454 False True True 1.05 1.96 False
roi bg both unpaired 3369 820 -2.58 t=2.14 p=0.0643 4570 10596 t=2.99 p=0.0173 False True True 1.05 1.86 False
roi bg both paired 3369 820 -2.58 t=2.02 p=0.114 4570 10596 t=2.81 p=0.0482 False True True 1.05 1.86 False
roi bg e1+e2 unpaired 3369 820 -2.58 t=2.14 p=0.0643 4570 10675 t=3.11 p=0.0145 False True True 1.05 1.97 False
roi bg e1+e2 paired 3369 820 -2.58 t=2.02 p=0.114 4570 10675 t=2.87 p=0.0454 False True True 1.05 1.97 False
loess0.1 raw both unpaired 2778 1134 -2.44 t=3.32 p=0.0105 5527 10910 t=2.67 p=0.0283 True True False 2.19 1.55 True
loess0.1 raw both paired 2778 1134 -2.44 t=3.21 p=0.0327 5527 10910 t=2.51 p=0.066 True False False 2.19 1.55 False
loess0.1 raw e1+e2 unpaired 2778 1134 -2.44 t=3.32 p=0.0105 5527 10989 t=2.78 p=0.0239 True True False 2.19 1.65 True
loess0.1 raw e1+e2 paired 2778 1134 -2.44 t=3.21 p=0.0327 5527 10989 t=2.63 p=0.0582 True False False 2.19 1.65 False
loess0.1 bg both unpaired 2566 922 -2.44 t=3.34 p=0.0102 5213 10596 t=2.68 p=0.0278 True True False 2.21 1.56 True
loess0.1 bg both paired 2566 922 -2.44 t=3.21 p=0.0327 5213 10596 t=2.51 p=0.066 True False False 2.21 1.56 False
loess0.1 bg e1+e2 unpaired 2566 922 -2.44 t=3.34 p=0.0102 5213 10675 t=2.79 p=0.0234 True True False 2.21 1.67 True
loess0.1 bg e1+e2 paired 2566 922 -2.44 t=3.21 p=0.0327 5213 10675 t=2.63 p=0.0582 True False False 2.21 1.67 False
loess0.2 raw both unpaired 2522 1181 -2.65 t=3.74 p=0.00567 5736 10910 t=2.56 p=0.0338 True True False 2.62 1.41 True
loess0.2 raw both paired 2522 1181 -2.65 t=3.34 p=0.0289 5736 10910 t=2.25 p=0.0879 True False False 2.62 1.41 False
loess0.2 raw e1+e2 unpaired 2522 1181 -2.65 t=3.74 p=0.00567 5736 10989 t=2.66 p=0.0287 True True False 2.62 1.52 True
loess0.2 raw e1+e2 paired 2522 1181 -2.65 t=3.34 p=0.0289 5736 10989 t=2.35 p=0.0785 True False False 2.62 1.52 False
loess0.2 bg both unpaired 2310 969 -2.65 t=3.77 p=0.00544 5422 10596 t=2.56 p=0.0334 True True False 2.64 1.42 True
loess0.2 bg both paired 2310 969 -2.65 t=3.34 p=0.0289 5422 10596 t=2.25 p=0.0879 True False False 2.64 1.42 False
loess0.2 bg e1+e2 unpaired 2310 969 -2.65 t=3.77 p=0.00544 5422 10675 t=2.67 p=0.0283 True True False 2.64 1.53 True
loess0.2 bg e1+e2 paired 2310 969 -2.65 t=3.34 p=0.0289 5422 10675 t=2.35 p=0.0785 True False False 2.64 1.53 False
loess0.3 raw both unpaired 2369 1198 -2.71 t=3.81 p=0.00519 5781 10910 t=2.54 p=0.0347 True True False 2.70 1.41 True
loess0.3 raw both paired 2369 1198 -2.71 t=3.32 p=0.0294 5781 10910 t=2.26 p=0.0862 True False False 2.70 1.41 False
loess0.3 raw e1+e2 unpaired 2369 1198 -2.71 t=3.81 p=0.00519 5781 10989 t=2.65 p=0.0294 True True False 2.70 1.51 True
loess0.3 raw e1+e2 paired 2369 1198 -2.71 t=3.32 p=0.0294 5781 10989 t=2.37 p=0.0769 True False False 2.70 1.51 False
loess0.3 bg both unpaired 2157 986 -2.71 t=3.84 p=0.00494 5467 10596 t=2.55 p=0.0342 True True False 2.73 1.41 True
loess0.3 bg both paired 2157 986 -2.71 t=3.32 p=0.0294 5467 10596 t=2.26 p=0.0862 True False False 2.73 1.41 False
loess0.3 bg e1+e2 unpaired 2157 986 -2.71 t=3.84 p=0.00494 5467 10675 t=2.66 p=0.029 True True False 2.73 1.52 True
loess0.3 bg e1+e2 paired 2157 986 -2.71 t=3.32 p=0.0294 5467 10675 t=2.37 p=0.0769 True False False 2.73 1.52 False
1P x 40 centre > trough holds in 24 of 32 readings; 1T x 40 periphery > centre in 20 ; both in 12
1T x 40 centre is the lowest point of its mean profile (reported, not required): {('roi', True): 8, ('loess0.1', False): 8, ('loess0.2', False): 8, ('loess0.3', False): 8}
1T x 40, O4 raw: centre 4884 to 5781 a.u., most peripheral value 10910 to 10989 a.u. (main geometry)
1T x 40, O4 bg: centre 4570 to 5467 a.u., most peripheral value 10596 to 10675 a.u. (main geometry)
1T x 40: mean profile minimum at -0.00 mm (4884 a.u.), centre 4884
All readings under every ROI geometry, split by paired and unpaired test (O6):
end width trough test readings holds_P holds_T holds max_p_P max_p_T directions_as_claimed footprint 1 group mean paired 16 12 4 0 0.11400 0.0879 True footprint 1 group mean unpaired 16 12 16 12 0.06470 0.0347 True footprint 1 leave one out paired 16 12 4 0 0.11600 0.0879 True footprint 1 leave one out unpaired 16 12 16 12 0.09040 0.0347 True footprint 3 group mean paired 16 16 2 2 0.03440 0.0804 True footprint 3 group mean unpaired 16 16 16 16 0.01290 0.0305 True footprint 3 leave one out paired 16 12 2 0 0.05870 0.0804 True footprint 3 leave one out unpaired 16 16 16 16 0.02160 0.0305 True footprint 5 group mean paired 16 16 0 0 0.03690 0.0904 True footprint 5 group mean unpaired 16 16 16 16 0.00727 0.0329 True footprint 5 leave one out paired 16 16 0 0 0.03840 0.0904 True footprint 5 leave one out unpaired 16 16 16 16 0.01370 0.0329 True needle 1 group mean paired 16 12 4 0 0.11600 0.0877 True needle 1 group mean unpaired 16 12 16 12 0.06930 0.0292 True needle 1 leave one out paired 16 12 4 0 0.14000 0.0877 True needle 1 leave one out unpaired 16 12 16 12 0.09510 0.0292 True needle 3 group mean paired 16 16 10 10 0.03450 0.0655 True needle 3 group mean unpaired 16 16 16 16 0.01480 0.0191 True needle 3 leave one out paired 16 12 10 6 0.05950 0.0655 True needle 3 leave one out unpaired 16 16 16 16 0.02370 0.0191 True needle 5 group mean paired 16 16 10 10 0.02630 0.0691 True needle 5 group mean unpaired 16 16 16 16 0.00760 0.0204 True needle 5 leave one out paired 16 16 10 10 0.04180 0.0691 True needle 5 leave one out unpaired 16 16 16 16 0.01190 0.0204 True Geometries where both comparisons hold under every reading (paired and unpaired): [] Per dish (raw ROI, averaged diagonals): centre, mean of both ends, at widths 1, 3 and 5 (measured footprint) 1P x 40 dish15 1827 2498 2038 2161 1959 2244 1P x 40 dish2 7941 1575 4912 1442 4240 1410 1P x 40 dish2-2 4280 1178 3235 1112 3112 1130 1P x 40 dish4 1985 1796 2157 2050 2129 2084 1P x 40 dish6 1872 1236 1997 1225 2147 1277 1T x 40 dish19 4776 15827 6098 14720 6039 13902 1T x 40 dish3 4714 6264 5249 6609 4908 6916 1T x 40 dish5 4902 14825 4871 14374 4716 14743 1T x 40 dish5-2 6136 6664 5964 6618 6280 6297 1T x 40 dish7 3891 10970 5254 11093 5789 10897
Result: not reproduced. No ROI geometry gives both comparisons under every reading. In the main geometry the 1P x 40 centre exceeds its trough in 24 of 32 readings, all LOESS readings (p = 0.00494 to 0.0327), not with the centre ROI (3581 against 1032 a.u., p = 0.0647 unpaired, 0.114 paired); the 1T x 40 centre is below its periphery in 20 of 32 (centre 4884 to 5781 a.u. raw and 4570 to 5467 after background, against 10910 to 10989 and 10596 to 10675), every unpaired test (p = 0.0145 to 0.0347) and the paired centre-ROI ones (p = 0.0454 and 0.0482), not the paired LOESS ones (p = 0.0582 to 0.0879). Both hold in 12 of 32. With ROIs 3 or 5 steps wide, both comparisons hold under all 16 unpaired readings (largest p 0.0329 for 1T x 40, 0.0148 for 1P x 40), but the paired 1T x 40 test holds in 0 to 10 of 16 (largest p 0.0655 to 0.0904). In every reading and geometry both differences go the way the paper says. The criterion asks for both comparisons "under all readings of O1, O2 and O4" with "a two-sided t-test (O6)", and O6 runs the paired and the unpaired test for comparisons inside the same dishes; this attempt reads the criterion as requiring both tests, as the fig5a-packet-center-enhanced criterion says outright. Under the other reading, in which the unpaired test alone counts, the outcome is the same with ROIs one step wide (both hold in 12 of 16 unpaired readings; the 1P x 40 centre ROI fails, p = 0.0643 to 0.0693 under either ROI end), and both comparisons hold with ROIs 3 or 5 steps wide. One step is the most direct reading of the ROI layout (see the section on the ROI width below), so under it the outcome does not depend on the reading of O6. The 1P x 40 centre-ROI mean rests mostly on dish2 (centre 7941 against ends 1575 a.u.; the other four dishes 1827 to 4280); in its full-resolution crop (section ROI geometry below) the centre ROI has a mean of 8151 a.u. with no saturated pixel, and pixels above 3 times the crop median cover 44.6% of it in two connected regions, the larger 49 um across, against 1.4 to 8.2% in the other four dishes. By eye the regions are clusters of round spots of cell size, consistent with a patch of labelled cells; the measure does not rule out debris. The 1T x 40 centre is the lowest point of its ROI mean profile (4884 a.u.), not of its LOESS fits.
fig5c-normalized-contrast¶
Normalized to the most peripheral ROI (100%), rotated Packet 2 and rotated trains show sharply different spatial distributions. The packet's center is relatively high and the trains' center is relatively low.
Criterion (claims file): the normalized centre (O1, O2, O3, O4) of 1P x 40 above that of 1T x 40 by an unpaired two-sided t-test at p < 0.05, under all readings; the normalized values, and whether 1P x 40 is above 100%, are reported. O3 (a): each dish normalized to its own peripheral ROI, then averaged; (b): the group's mean centre over its mean periphery. For the test under O3 (b) each dish's centre is divided by its group's mean periphery: an assumption of this attempt, as the claims file gives the reported value but not the per-dish values a test needs.
norm_profiles = {}
R5c = test_5c(keep=norm_profiles)
with pd.option_context("display.width", 250, "display.max_columns", 20):
print(R5c.to_string(index=False))
print("\nReadings where 1P x 40 > 1T x 40 at p < 0.05:", int(R5c.holds.sum()), "of", len(R5c),
"; 1P x 40 normalized centre above 100% in", int(R5c.P_above_100.sum()))
print("Ranges: 1P x 40 %.0f to %.0f %%, 1T x 40 %.0f to %.0f %%; p from %.2g to %.2g" % (
R5c.P_norm_centre.min(), R5c.P_norm_centre.max(), R5c.T_norm_centre.min(), R5c.T_norm_centre.max(), R5c.p.min(), R5c.p.max()))
print("\nAll readings under every ROI geometry:")
rows = []
for g in GEOMS:
R = test_5c(g)
fail = R[~R.holds]
rows.append(dict(end=g[0], width=g[1], readings=len(R), holds=int(R.holds.sum()), P_above_100=int(R.P_above_100.sum()),
P_range=f"{R.P_norm_centre.min():.0f}-{R.P_norm_centre.max():.0f}", T_range=f"{R.T_norm_centre.min():.0f}-{R.T_norm_centre.max():.0f}",
max_p=float(f"{R.p.max():.3g}"), failing=", ".join(f"{a}/{b}/{c}/{e} p={p:.3g}" for a, b, c, e, p in
fail[["O1", "O4", "O2", "O3", "p"]].values)))
G5c = pd.DataFrame(rows)
with pd.option_context("display.width", 300, "display.max_columns", 20, "display.max_colwidth", 200):
print(G5c.to_string(index=False))
fig, ax = plt.subplots(1, 2, figsize=(13, 4))
for a, o4 in zip(ax, ["raw", "bg"]):
for name, c in [("P", "C0"), ("T", "C3")]:
Y = norm_profiles[(o4, name)]
a.errorbar(x_mm, Y.mean(0), Y.std(0, ddof=1) / np.sqrt(len(Y)), color=c, lw=1,
label=("1P x 40" if name == "P" else "1T x 40") + f" (range {Y.mean(0).min():.0f} to {Y.mean(0).max():.0f} %)")
a.axhline(100, color="k", lw=0.5); a.set(title=f"Fig. 5C, O4 {o4}, O3 (a), O2 both ends", xlabel="distance from centre (mm)", ylabel="% of most peripheral ROI")
a.legend(fontsize=8)
plt.tight_layout(); plt.show()
O1 O4 O2 O3 P_norm_centre T_norm_centre diff_points t p holds P_above_100 gap
roi raw both a 240.5 53.2 187.3 2.23 0.056000 False True 1.11
roi raw both b 216.2 44.8 171.4 2.40 0.043500 True True 1.35
roi raw e1+e2 a 220.9 53.0 167.9 2.27 0.052900 False True 1.13
roi raw e1+e2 b 193.1 44.4 148.7 2.33 0.048500 True True 1.28
roi bg both a 266.6 51.9 214.7 2.19 0.060100 False True 1.07
roi bg both b 233.2 43.1 190.1 2.32 0.048800 True True 1.28
roi bg e1+e2 a 242.3 51.7 190.6 2.23 0.056200 False True 1.10
roi bg e1+e2 b 205.1 42.8 162.3 2.25 0.054300 False True 1.21
loess0.1 raw both a 184.3 59.7 124.7 2.93 0.019100 True True 1.75
loess0.1 raw both b 167.7 50.7 117.0 3.98 0.004070 True True 2.90
loess0.1 raw e1+e2 a 170.9 59.1 111.8 2.80 0.023100 True True 1.62
loess0.1 raw e1+e2 b 149.8 50.3 99.5 3.78 0.005360 True True 2.70
loess0.1 bg both a 200.2 58.6 141.6 2.84 0.021700 True True 1.68
loess0.1 bg both b 177.6 49.2 128.4 3.83 0.005050 True True 2.75
loess0.1 bg e1+e2 a 184.1 58.0 126.0 2.73 0.025800 True True 1.55
loess0.1 bg e1+e2 b 156.2 48.8 107.4 3.63 0.006650 True True 2.55
loess0.2 raw both a 167.1 63.3 103.8 2.84 0.021800 True True 1.67
loess0.2 raw both b 152.2 52.6 99.6 4.81 0.001340 True True 3.71
loess0.2 raw e1+e2 a 155.4 62.7 92.8 2.59 0.032300 True True 1.39
loess0.2 raw e1+e2 b 136.0 52.2 83.8 4.51 0.001980 True True 3.40
loess0.2 bg both a 180.0 62.4 117.6 2.77 0.024400 True True 1.58
loess0.2 bg both b 159.9 51.2 108.7 4.60 0.001750 True True 3.50
loess0.2 bg e1+e2 a 166.1 61.8 104.4 2.52 0.035900 True True 1.31
loess0.2 bg e1+e2 b 140.6 50.8 89.9 4.31 0.002580 True True 3.20
loess0.3 raw both a 157.2 63.5 93.7 2.74 0.025300 True True 1.56
loess0.3 raw both b 143.0 53.0 90.0 5.13 0.000891 True True 4.04
loess0.3 raw e1+e2 a 146.5 62.9 83.6 2.44 0.040300 True True 1.23
loess0.3 raw e1+e2 b 127.8 52.6 75.1 4.78 0.001390 True True 3.67
loess0.3 bg both a 168.4 62.6 105.8 2.66 0.028600 True True 1.47
loess0.3 bg both b 149.3 51.6 97.7 4.89 0.001210 True True 3.79
loess0.3 bg e1+e2 a 155.9 62.0 93.9 2.37 0.045100 True True 1.15
loess0.3 bg e1+e2 b 131.3 51.2 80.1 4.54 0.001890 True True 3.44
Readings where 1P x 40 > 1T x 40 at p < 0.05: 27 of 32 ; 1P x 40 normalized centre above 100% in 32
Ranges: 1P x 40 128 to 267 %, 1T x 40 43 to 64 %; p from 0.00089 to 0.06
All readings under every ROI geometry:
end width readings holds P_above_100 P_range T_range max_p failing footprint 1 32 27 32 128-267 43-64 0.0601 roi/raw/both/a p=0.056, roi/raw/e1+e2/a p=0.0529, roi/bg/both/a p=0.0601, roi/bg/e1+e2/a p=0.0562, roi/bg/e1+e2/b p=0.0543 footprint 3 32 32 32 132-218 50-66 0.0376 footprint 5 32 31 32 129-204 51-67 0.0505 loess0.3/bg/e1+e2/a p=0.0505 needle 1 32 27 32 164-352 46-69 0.0824 roi/raw/both/a p=0.0619, roi/raw/e1+e2/a p=0.0752, roi/bg/both/a p=0.0696, roi/bg/e1+e2/a p=0.0824, roi/bg/e1+e2/b p=0.0512 needle 3 32 32 32 145-249 52-64 0.0489 needle 5 32 32 32 134-218 52-63 0.0331
Result: not reproduced, resting on ROIs one step wide (assumed). The 1P x 40 normalized centre is above 100% under all 32 readings in every geometry (128 to 267% in the main one) and the 1T x 40 one below (43 to 64%), so the direction holds everywhere. Whether the difference is significant under every reading, as the criterion requires, depends on the ROI width, which the paper does not give. ROIs one step wide, the square ROIs that tile the diagonal at the spacing of Methods 2.4, are the most direct reading of the paper's layout (next section), and this attempt takes that width. With it the difference is significant in 27 of 32 readings, under either reading of the ROI end; the 5 that fail use the centre ROI (p = 0.0512 to 0.0824), so the claim is not reproduced, a near miss. With ROIs three steps wide it holds in all 32 under both ROI ends (largest p 0.0376 and 0.0489), and with five steps in 31 of 32 (one at p = 0.0505) or all 32 (largest p 0.0331). The claim depends on fig5b-rotation-center-peak, which is not reproduced.
The ROI width of Methods 2.4¶
Methods 2.4 says the profiles used "92 or 224 regions of interest (ROI), placed along one or both diagonals", over +/- 3.65 mm from the centre. 224 ROIs on two diagonals are 112 per 7.3-mm diagonal, one every 65 um. ROIs that follow each other along the diagonal without gaps or overlap at that spacing are about 65 um long; the most direct reading of that layout is a square ROI, as wide as it is long: the width of one step. Below, fig5b-rotation-center-peak and fig5c-normalized-contrast under that width, with both readings of the ROI end, and for fig5b with the unpaired test alone as well (the reading of its criterion in which the paired test of O6 is not required).
steps = np.array([1000 * d["mm"][:, 1].mean() for d in D])
print("Methods 2.4: 112 ROIs per 2 x 3.65-mm diagonal, one every %.1f um; ROI step measured here %.1f to %.1f um (at 0.65 um per pixel)"
% (2 * 3.65e3 / 112, steps.min(), steps.max()))
rows = []
for end in cancan.ROI_ENDS:
b, c = test_5b((end, 1)), test_5c((end, 1))
u = b[b.test == "unpaired"]
fail_u = u[~u.holds]
rows.append(dict(end=end, width=1, fig5b_all=f"{int(b.holds.sum())} of {len(b)}", fig5b_unpaired=f"{int(u.holds.sum())} of {len(u)}",
fig5b_unpaired_failing=", ".join(sorted({f"{a}/{o4}: 1P p={p:.3g}" if not hp else f"{a}/{o4}/{o2}: 1T p={pt:.3g}"
for a, o4, o2, p, pt, hp in fail_u[["O1", "O4", "O2", "p_P", "p_T", "holds_P"]].values})),
fig5c=f"{int(c.holds.sum())} of {len(c)}", fig5c_failing_p=f"{c[~c.holds].p.min():.3g} to {c[~c.holds].p.max():.3g}" if (~c.holds).any() else "-"))
with pd.option_context("display.width", 300, "display.max_columns", 20, "display.max_colwidth", 200):
print(pd.DataFrame(rows).to_string(index=False))
Methods 2.4: 112 ROIs per 2 x 3.65-mm diagonal, one every 65.2 um; ROI step measured here 65.4 to 68.4 um (at 0.65 um per pixel)
end width fig5b_all fig5b_unpaired fig5b_unpaired_failing fig5c fig5c_failing_p footprint 1 12 of 32 12 of 16 roi/bg: 1P p=0.0643, roi/raw: 1P p=0.0647 27 of 32 0.0529 to 0.0601 needle 1 12 of 32 12 of 16 roi/bg: 1P p=0.0688, roi/raw: 1P p=0.0693 27 of 32 0.0512 to 0.0824
ROI width. With ROIs one step wide, fig5c-normalized-contrast fails in 5 of 32 readings under either ROI end, all of them with the centre ROI (p = 0.0529 to 0.0601 with the measured footprint, 0.0512 to 0.0824 with the needle reading). After round 2 this attempt proposed not_reproduced for it, resting on that width (assumed), in place of the verified blocked (methods_underspecified); the third verification round agreed, and that is now its outcome. Under the same width fig5b-rotation-center-peak fails whether or not the paired test is required: with the unpaired test alone both comparisons hold in 12 of 16 readings, and the 4 that fail are the 1P x 40 centre ROI (p = 0.0643 to 0.0693). Its verified outcome, not reproduced, stands under either reading of its criterion at this width, so no change was proposed for it.
ROI geometry and the 1P x 40 centre¶
The paper gives no ROI width and places the ROIs by eye. Every claim above was run under six ROI geometries. The ROI end under the measured footprint reading moves with the footprint radius estimated from the image; the needle reading fixes it from the paper's geometry. Here: the ROI ends of both readings per group, and, because the 1P x 40 centre ROI mean is dominated by one dish (dish2), a full-resolution crop of the array centre of each 1P x 40 dish, to see whether that value comes from cells or from debris.
from scipy import ndimage as ndi
rows = []
for g in groups:
for d in group(g):
px_mm = cancan.BLOCK * cancan.UM_PER_PIXEL / 1000
rows.append(dict(group=g, dish=d["dish"], footprint_radius_mm=d["el"]["radius"] * px_mm,
end_footprint_mm=d["mm"][:, -1].mean(),
end_needle_mm=np.mean([cancan.roi_positions(d["el"], i, "needle")[1] for i in range(4)]) * cancan.N_HALF * px_mm))
E = pd.DataFrame(rows)
print(E.groupby("group")[["footprint_radius_mm", "end_footprint_mm", "end_needle_mm"]].agg(["min", "max"]).round(3).to_string())
print("Largest shift of the ROI end between readings: %.3f mm (%.1f ROI steps)" % (
(E.end_footprint_mm - E.end_needle_mm).abs().max(), (E.end_footprint_mm - E.end_needle_mm).abs().max() / (E.end_footprint_mm.mean() / 56)))
HALF = 400 # full-resolution pixels, 0.26 mm at 0.65 um per pixel
step_px = d["mm"][:, 1].mean() * 1000 / cancan.UM_PER_PIXEL # one ROI step in full-resolution pixels
fig, axs = plt.subplots(1, len(RP), figsize=(4.2 * len(RP), 4.4))
for a, d in zip(axs, RP):
c = d["el"]["centre"] * cancan.BLOCK + (cancan.BLOCK - 1) / 2
x = cancan.load_crop(d["url"], c, HALF).astype(float)
h = int(step_px / 2)
roi = x[HALF - h:HALF + h + 1, HALF - h:HALF + h + 1]
# bright objects in the centre ROI: connected regions above 3 times the crop median, of at least 20 pixels
lab, n = ndi.label(roi > 3 * np.median(x))
area = np.bincount(lab.ravel())[1:]
area = area[area >= 20]
diam = 2 * np.sqrt(area / np.pi) * cancan.UM_PER_PIXEL
print(f"1P x 40 {d['dish']:8s}: crop {x.shape}, crop mean {x.mean():.0f}, median {np.median(x):.0f}, "
f"centre ROI square ({2 * h + 1} px) mean {roi.mean():.0f}, its 99th percentile {np.percentile(roi, 99):.0f}, "
f"pixels >= 65000: {int((x >= 65000).sum())}; bright objects in the ROI: {len(area)}, "
f"equivalent diameter {np.median(diam) if len(diam) else np.nan:.1f} um (median), largest {diam.max() if len(diam) else np.nan:.1f} um, "
f"{100 * area.sum() / roi.size:.1f}% of the ROI")
a.imshow(x, cmap="gray", vmin=0, vmax=np.percentile(x, 99.8))
a.add_patch(plt.Rectangle((HALF - h, HALF - h), 2 * h, 2 * h, fill=False, color="C1"))
a.set_title(f"1P x 40 {d['dish']}", fontsize=10); a.axis("off")
plt.tight_layout(); plt.show()
print("Each panel is 0.52 mm across at 0.65 um per pixel; brightness to the 99.8th percentile of the crop. Orange: a square one ROI step "
"wide on the array centre (the analysis ROI has the same size, turned along the diagonal).")
footprint_radius_mm end_footprint_mm end_needle_mm
min max min max min max
group
1P x 40 0.575 0.660 3.660 3.829 3.583 3.672
1T x 40 0.627 0.720 3.737 3.799 3.665 3.714
40 Packets 0.627 0.666 3.666 3.740 3.594 3.641
40 Trains 0.634 0.702 3.677 3.762 3.604 3.679
Largest shift of the ROI end between readings: 0.156 mm (2.3 ROI steps)
1P x 40 dish15 : crop (801, 801), crop mean 2202, median 1711, centre ROI square (105 px) mean 1861, its 99th percentile 5644, pixels >= 65000: 671; bright objects in the ROI: 4, equivalent diameter 4.4 um (median), largest 5.7 um, 1.4% of the ROI
1P x 40 dish2 : crop (801, 801), crop mean 3632, median 2531, centre ROI square (105 px) mean 8151, its 99th percentile 25293, pixels >= 65000: 0; bright objects in the ROI: 2, equivalent diameter 32.3 um (median), largest 49.0 um, 44.6% of the ROI
1P x 40 dish2-2 : crop (801, 801), crop mean 2858, median 2377, centre ROI square (105 px) mean 3816, its 99th percentile 29206, pixels >= 65000: 136; bright objects in the ROI: 4, equivalent diameter 9.3 um (median), largest 16.6 um, 8.2% of the ROI
1P x 40 dish4 : crop (801, 801), crop mean 1978, median 1646, centre ROI square (105 px) mean 1941, its 99th percentile 5814, pixels >= 65000: 0; bright objects in the ROI: 5, equivalent diameter 6.1 um (median), largest 6.4 um, 2.6% of the ROI
1P x 40 dish6 : crop (801, 801), crop mean 1949, median 1625, centre ROI square (105 px) mean 1940, its 99th percentile 12049, pixels >= 65000: 0; bright objects in the ROI: 3, equivalent diameter 7.1 um (median), largest 10.6 um, 3.2% of the ROI
Each panel is 0.52 mm across at 0.65 um per pixel; brightness to the 99.8th percentile of the crop. Orange: a square one ROI step wide on the array centre (the analysis ROI has the same size, turned along the diagonal).
Helper checks on a synthetic image¶
The helpers run on a synthetic image of the same size as the block images, with four dark footprints of known centre and radius on a linear intensity ramp, where the mean of any ROI equals the ramp at its centre. Also the profile indexing, the LOESS, the t-tests and the error-bar rule against values computed by hand.
# Helper checks on a synthetic image with known geometry and intensity
rng = np.random.default_rng(0)
N, rad = 2790, 100
yy, xx = np.mgrid[:N, :N].astype(np.float32)
img = 1000 + 2 * yy # linear ramp: a band mean equals the value at its centre
true = np.array([[723, 723], [723, 2067], [2067, 2067], [2067, 723]], float) # clockwise from top left
for r0, c0 in true:
img[(yy - r0) ** 2 + (xx - c0) ** 2 <= rad ** 2] = 200 # dark footprints
img += rng.normal(0, 5, img.shape).astype(np.float32)
el = cancan.find_electrodes(img)
print("footprint centres found (row, col):", el["centres"].round(1).tolist(), "; largest error %.1f blocks (a disk of radius 85 fits anywhere within 15 blocks of the true centre)" % np.abs(el["centres"] - true).max())
print("footprint radius found %.1f (true %d); background %.1f (true 200)" % (el["radius"], rad, cancan.footprint_background(img, el)))
for end in cancan.ROI_ENDS:
for w in cancan.ROI_WIDTHS:
pr = cancan.diagonal_profiles(img, el, across=w, end=end)
pts = np.array([cancan.roi_positions(el, i, end)[0] for i in range(4)])
err = np.abs(pr["halves"] - (1000 + 2 * pts[..., 0])).max()
print(f" ROIs {end}, width {w}: largest |band mean - ramp at ROI centre| = {err:.2f} a.u.; outermost ROI at {pr['frac'][-1]:.2f} of the range")
y = np.arange(113, dtype=float)
h = np.array([np.arange(57.0) + 100 * i for i in range(4)]); h[:, 0] = 0 # the four halves share the centre ROI
sp = cancan.signed_profile(h, 0, 2)
print("signed_profile: index 0, 56, 112 =", sp[[0, 56, 112]], "(expected 256 = electrode 2 end, 0 = centre, 56 = electrode 0 end)")
print("LOESS of a line returns the line:", all(np.allclose(cancan.loess(3 * y + 7, s), 3 * y + 7) for s in cancan.LOESS_SPANS),
"; centre values", [round(cancan.centre_value(3 * y + 7, o), 3) for o in ["roi", "loess0.1", "loess0.3"]])
a, b = np.array([1.0, 2, 3, 4]), np.array([2.0, 4, 5, 9])
sp_ = np.sqrt(((a - a.mean()) ** 2).sum() / 6 + ((b - b.mean()) ** 2).sum() / 6)
t_manual = (a.mean() - b.mean()) / (sp_ * np.sqrt(0.5))
dd = a - b
print("ttest unpaired t = %.4f (by hand %.4f), paired t = %.4f (by hand %.4f)" % (
cancan.ttest(a, b)[0], t_manual, cancan.ttest(a, b, paired=True)[0], dd.mean() / (dd.std(ddof=1) / 2)))
sa, sb = a.std(ddof=1) / 2, b.std(ddof=1) / 2
print("errorbar_gap = %.4f (by hand %.4f)" % (cancan.errorbar_gap(a, b)[0], (abs(a.mean() - b.mean()) - sa - sb) / max(sa, sb)))
footprint centres found (row, col): [[735.0, 714.0], [727.0, 2079.0], [2078.0, 2057.0], [2053.0, 721.0]] ; largest error 14.0 blocks (a disk of radius 85 fits anywhere within 15 blocks of the true centre) footprint radius found 99.0 (true 100); background 200.0 (true 200) ROIs footprint, width 1: largest |band mean - ramp at ROI centre| = 1.64 a.u.; outermost ROI at 1.00 of the range ROIs footprint, width 3: largest |band mean - ramp at ROI centre| = 1.03 a.u.; outermost ROI at 1.00 of the range ROIs footprint, width 5: largest |band mean - ramp at ROI centre| = 0.57 a.u.; outermost ROI at 1.00 of the range ROIs needle, width 1: largest |band mean - ramp at ROI centre| = 1.93 a.u.; outermost ROI at 1.00 of the range ROIs needle, width 3: largest |band mean - ramp at ROI centre| = 1.13 a.u.; outermost ROI at 1.00 of the range ROIs needle, width 5: largest |band mean - ramp at ROI centre| = 0.55 a.u.; outermost ROI at 1.00 of the range signed_profile: index 0, 56, 112 = [256. 0. 56.] (expected 256 = electrode 2 end, 0 = centre, 56 = electrode 0 end) LOESS of a line returns the line: True ; centre values [175.0, 175.0, 175.0] ttest unpaired t = -1.5554 (by hand -1.5554), paired t = -2.8868 (by hand -2.8868) errorbar_gap = 0.2599 (by hand 0.2599)
Claims not tested¶
- fig3-field-simulation (Fig. 3, blocked,
data_missing). The values are outputs of the paper's Sim4Life model. Neither this dandiset nor any other dandiset citing the paper holds field maps, and the paper names no deposit. The criterion forbids a model built here. - fig4-sham-no-uptake, fig4-unipolar-uptake, fig4-bipolar-weaker, fig4-8phase-cancellation
(Fig. 4, blocked,
data_missing). They need the confocal time-lapse of single cells on coverslips (Methods 2.2). This dandiset holds only the monolayer images of Fig. 5, and none of the four other dandisets that cite the paper (000875, 001025, 001027, 001028) holds single-cell time-lapse data; a search of DANDI for the lab's other dandisets found none that cites the paper. fig4-8phase-cancellation depends on the other two Fig. 4 claims. - Figs. 6 to 9 (eleven claims, blocked,
data_missing). The data are not in dandiset 000686, the dataset the claims file names for them (data: dataset), so outside data may not be used (protocol, author, step 4). Four other dandisets of the paper's authors, DANDI:000875 (Fig. 6), 001025 (Fig. 7), 001027 (Fig. 8) and 001028 (Fig. 9), each name this paper and appear to hold them; they were not opened here. A question asks a reader whether these claims should bedata: elsewhere(questions/000686.yaml).
Summary¶
| Claim | Outcome |
|---|---|
| fig5a-trains-stronger | not reproduced: trains do not exceed packets at the e3 end under any reading or ROI geometry |
| fig5a-packet-center-enhanced | not reproduced: e1 end above centre in all 4 dishes, not significant under any reading or ROI geometry |
| fig5b-rotation-center-peak | not reproduced: directions as claimed, but no ROI geometry gives significance under every reading (paired 1T x 40 test p up to 0.09); with ROIs one step wide the 1P x 40 centre ROI fails with the unpaired test too |
| fig5c-normalized-contrast | not reproduced, assuming ROIs 1 step wide (the most direct reading of Methods 2.4): direction as claimed in every reading, but the centre ROI misses p < 0.05 in 5 of 32 readings (p = 0.051 to 0.082); with ROIs 3 steps wide it holds under every reading |
| fig3-field-simulation | blocked, data_missing |
| four Fig. 4 claims | blocked, data_missing |
| eleven Fig. 6 to 9 claims | blocked, data_missing: not in this dandiset (apparently in DANDI:000875, 001025, 001027, 001028) |
A failure to reproduce here may be a failure of this attempt, not of the paper: see the reasoning and the differences from the paper's method.
How the analysis differs from the paper's. The paper's ROI measurements were made in MetaMorph and are not deposited; here the electrodes are found automatically and ROIs are placed by the layout of Methods 2.4, with ROIs one step long along the diagonal and 1, 3 or 5 steps wide, ending at the measured footprint or at a fixed fraction of the electrode distance (the paper gives neither ROI size nor how the ROIs were placed). With n = 4 or 5 dishes the tests on the single centre ROI are close to p = 0.05 and change with the ROI width; this decides fig5c-normalized-contrast, not the other three, and with ROIs one step wide, the most direct reading of the layout, it fails. e1 is not labelled in the files; the electrode brightest after non-rotated treatment is the top-right one in all 8 dishes, which is also the mirror image of Fig. 3's top-left electrode. The background of reading O4 (b) is the median inside the electrode footprints; it changes no outcome. Distances in mm assume 0.65 um per pixel, which the files do not state.
What was found about the dataset. It holds only Fig. 5's 18 dishes; the other monolayer experiments appear, from their titles and descriptions, to be in four other dandisets that cite the paper (not opened here). The species is given as Bos taurus (NCBITaxon_9913), where the paper used CHO-K1 hamster cells; the subject age is P1D for a cell line; the dandiset's assets summary gives 36 subjects for 18 files with 18 distinct subject ids; the image has no pixel size (resolution -1) and no unit; the electrodes are not labelled; the file description speaks of a six-well plate where the paper used 35-mm dishes. Every subject id is unique and no two images are near duplicates (highest correlation 0.815), so the files with the same dish number ("dish2" and "dish2-2", "dish11" in two groups) are different dishes, consistent with their session dates.