Reproducing Semenov, Kim, Silkuniene & Pakhomov (2025) from the DANDI datasets it names¶
Paper. Semenov I, Kim V, Silkuniene G, Pakhomov AG (2025). Excitation and polarization of isolated neurons by high-frequency sine waves for temporal interference stimulation. Cell Reports Physical Science 6, 102660. doi:10.1016/j.xcrp.2025.102660.
Text the claims file was written from. The NIH author manuscript (NIHMS2098259, PMC12320893) as served by Europe PMC, with its supplement (Figures S1 to S4, Table S1); the version of record was not read. The paper reports no original code.
Dataset. The target is DANDI:001048, version 0.250701.1624, its only published version, pinned
for that reason. It names the paper as related, but it holds TIRF images of HEK cells under patch
clamp from another study ("Plasma membrane conductance changes in response to supraphysiological
hyperpolarization"), shown below. The paper's own data are in the ten dandisets its Data and code
availability lists, and the claims file marks every claim data: elsewhere in them. Each is read at
its only published version, the most recent:
| Dandiset | Version | Holds | Used for |
|---|---|---|---|
| 001030 | 0.240520.1924 | 2-s ramped 2-kHz waves in NaGlu, 4 neurons | Figure 3C |
| 001031 | 0.240520.1924 | 100-kHz bipolar pulses | not in the paper, not used |
| 001032 | 0.240527.1432 | 250-ms waves, 14 neurons | Figures 1, 2 |
| 001033 | 0.240527.1432 | 250-ms waves, 10 neurons | Figures 1, 2 |
| 001034 | 0.240527.1431 | a copy of 001032 | duplicate check only |
| 001035 | 0.240527.1432 | 250-ms waves, 15 neurons | Figures 1, 2 |
| 001271 | 0.250805.1352 | TTX and vehicle | Figure S2 |
| 001272 | 0.250805.1358 | Di-8-ANEPPS, 3 neurons | Figure S1 |
| 001273 | 0.250805.1356 | strobe imaging of 2-us pulses | Figure 4 |
| 001274 | 0.250805.1401 | strobe imaging of pulses and sine waves | Figures 5, 6, S3 |
Every file is streamed from the archive; nothing is downloaded to disk. The paper does not say which figure is in which dandiset; the assignment is the claims file's, from the titles, and the data checks below are consistent with it.
Claims.
| Claim | Tested or blocked |
|---|---|
| ap-selection | tested |
| fig1-pure-sine-excites-all | tested |
| fig1-excitation-block | tested |
| fig1-20khz-needs-twice-the-field | computed; blocked, metadata insufficient (voltages not calibrated per frequency) |
| fig2d-5ap-above-1ap | tested |
| fig2d-modulated-thresholds-higher | tested |
| fig2-modulated-locks-to-beats | tested |
| fig2c-pure-more-aps | tested |
| fig3b-low-conductivity-halves-threshold | blocked, data missing (no NaGlu 250-ms sessions) |
| fig3-extrapolated-threshold | blocked, depends on fig3b |
| fig3-not-heating | computed; blocked, methods underspecified |
| fig3c-sustained-firing-2s | tested |
| figs1-di8-same-as-fluovolt | tested (3 of 8 neurons deposited) |
| figs2-ttx-blocks-depolarization | tested |
| fig4-pooled-time-constants | computed; blocked, metadata insufficient (no ROIs, experiments not identifiable) |
| fig4-orientation-and-poles | blocked, metadata insufficient (cell shape and orientation groups) |
| fig5-sine-asymmetry | blocked, metadata insufficient (Figure 5 and S3 cells not identifiable) |
| fig5-whole-cell-depolarization | blocked, metadata insufficient (same cells) |
| fig6-sine-less-than-pulse | computed; blocked, metadata insufficient (no ROIs, cells not identifiable, fields not measured) |
| threshold-depolarization-estimate | computed; blocked, metadata insufficient (voltages not calibrated) |
Revision after the first verification round. The verifier was given the claims file of another dandiset (001274) for the same paper, so its claim ids and readings differ from this file's; its points on the code apply here as well and were taken up. The AP detector now asks for a width of 2 to 20 ms both at half the peak's dF and at half its prominence (section Readings; the first round used the prominence only, which counted small oscillations on a raised plateau as APs), and the single measures are shown as variants at the end. Added: a check that APs within 5 ms of onset are not onset artefacts, the late dF of the modulated and TTX sessions, the APs that remain in the 2-s sessions, the carrier phase of the strobe oscillations and the pulse onset. The strobe LOESS is now R's loess. Three claims whose size depends on calibrations or ROIs the dandisets do not hold are now blocked.
Revision after the second verification round. Text corrected where the outputs did not support it. Added: the session-description wording on TTX and vehicle (Data check), the O4 levels of the modulated session that decides fig2-modulated-locks-to-beats, and a box-size variant of the strobe ROIs (fig6-sine-less-than-pulse). fig2d-modulated-thresholds-higher is now marked as resting on the crop ROI. For fig6-sine-less-than-pulse a change to blocked, metadata insufficient, is proposed for a third round, as for fig4-pooled-time-constants; its outcome is unchanged here.
Revision after the third verification round. The verifier agreed with the proposed change, and fig6-sine-less-than-pulse is now blocked, metadata insufficient. Added: the per-cycle reading with the full peak-to-peak amplitude, which the criterion names (the main table regresses half of it).
The notebook takes about 40 minutes, most of it streaming the strobe dandisets (about 11 GB of 001273 and the 001274 files of nine cells).
import time, warnings
warnings.filterwarnings("ignore", category=RuntimeWarning)
import numpy as np, pandas as pd
from scipy import stats
import matplotlib.pyplot as plt
import tisine as ts
pd.set_option("display.width", 200); pd.set_option("display.max_columns", 40); pd.set_option("display.max_rows", 200)
TL = ["001030", "001032", "001033", "001034", "001035", "001271", "001272"]
Load the data¶
Every time-lapse file is streamed and reduced, as it loads, to the soma dF traces under each ROI reading (O1 a, b and c) and F0 reading (O2 a and b). The strobe dandisets 001273 and 001274 are read later, file by file, in their own sections. Dandiset 001048, the one that names the paper, is listed and one of its files opened, to show what it holds.
m48 = ts.dandiset_meta("001048")
print("001048:", m48["name"]); print(" species:", [s["name"] for s in m48["assetsSummary"]["species"]])
print(" description:", m48["description"][:400])
a48 = ts.list_assets("001048")
for a in a48: print(" ", a["path"], round(a["size"] / 1e6), "MB")
import h5py, remfile
with h5py.File(remfile.File(ts.s3_url("001048", a48[0]["asset_id"])), "r") as f:
print(" one file:", list(f["acquisition"].keys()), f["acquisition/TIRFImaging/data"].shape,
"| strain:", f["general/subject/strain"][()].decode(), "| subject description:", f["general/subject/description"][()].decode())
print(" session description:", f["session_description"][()].decode()[:300])
001048: Plasma membrane conductance changes in response to supraphysiological hyperpolarization species: ['Homo sapiens - Human'] description: This study aimed to visualize and characterize membrane lesions caused by hyperpolarization and their protective roles against membrane overcharging. We identified three types of lesions in HEK cells, each with distinct characteristics. Hyperpolarization to specific voltages created diffuse zonal electropermeabilization, focal pores, and high-conductance pores, with adaptive conductance changes pr sub-Cell11/sub-Cell11_ses-session-11_image.nwb 674 MB sub-Cell03/sub-Cell03_ses-session-03_image.nwb 674 MB sub-Cell07/sub-Cell07_ses-session-07_image.nwb 688 MB sub-Cell12/sub-Cell12_ses-session-12_image.nwb 684 MB sub-Cell04/sub-Cell04_ses-session-04_image.nwb 751 MB sub-Cell16/sub-Cell16_ses-session-16_image.nwb 668 MB sub-Cell13/sub-Cell13_ses-session-13_image.nwb 650 MB
one file: ['TIRFImaging'] (1202, 800, 800) | strain: cellline: HEK 923 | subject description: Grown on Indium Tin Oxide (ITO) coated coverslip session description: Subject ID: Cell11. General description of the experimental setup and conditions. Bath solution was composed of (in mM): 136 tetraethylammonium chloride (TEACl), 2 MgCl2, 10 CaCl2, 10 HEPES, and 10 glucose (pH 7.3, 300-310 mOsm/kg). Pipette solution contained (in mM): 20 TEACl, 130 CsCl, 10 HEPES, 4
import time
t0 = time.time()
REC = {d: ts.build_records(d, invert=(d == "001272")) for d in TL}
for d in TL:
print(d, len(REC[d]), "files,", sum("error" in r for r in REC[d]), "failed to stream")
print("streamed in %.0f s" % (time.time() - t0))
TAB = {d: ts.session_table(REC[d], d) for d in TL}
001030 17 files, 0 failed to stream 001032 178 files, 0 failed to stream 001033 193 files, 0 failed to stream 001034 178 files, 0 failed to stream 001035 254 files, 0 failed to stream 001271 287 files, 0 failed to stream 001272 65 files, 0 failed to stream streamed in 468 s
Data check¶
# Contents of every dandiset the paper lists: titles, files, subjects, species, frame counts and rates.
for d in ["001030", "001031", "001032", "001033", "001034", "001035", "001271", "001272", "001273", "001274"]:
m = ts.dandiset_meta(d)
print(d, ts.VERSIONS[d], "|", m["name"][:100], "|", m["assetsSummary"]["numberOfFiles"], "files,",
round(m["assetsSummary"]["numberOfBytes"] / 1e9, 2), "GB, species", [s["name"] for s in m["assetsSummary"].get("species", [])])
for d in TL:
R = [r for r in REC[d] if "error" not in r]
print(f"\n{d}: {len(R)} files, {len(set(r['subject_id'] for r in R))} subject ids, species {sorted(set(r['species'] for r in R))}")
print(" subject descriptions:", sorted(set(r["subject_description"] for r in R)))
print(" frames per file:", pd.Series([r['n_frames'] for r in R]).value_counts().to_dict(), "| rates (frames/s):", sorted(set(r['rate'] for r in R)))
st = pd.DataFrame([r["stim"] for r in R])
print(" stimuli:", st.groupby(["kind", "carrier", "mod"], dropna=False).size().to_dict())
print(" TTX sessions:", int(st.ttx.sum()), "| labelled volts from", st.volts.min(), "to", st.volts.max())
print(" solution NaCl:", sum("140 NaCl" in r["session_description"] for r in R), "NaGluconate:", sum("NaGluconate" in r["session_description"] for r in R))
001030 0.240520.1924 | High-frequency (2 kHz) sine wave stimulation of individual hippocampal neurons: the effect of elect | 17 files, 0.05 GB, species ['Mus musculus - House mouse'] 001031 0.240520.1924 | High-frequency (100 kHz) rectangular bipolar pulses stimulation of individual hippocampal neurons: t | 10 files, 0.0 GB, species ['Mus musculus - House mouse'] 001032 0.240527.1432 | Characterization of Neuron Action Potential Triggering Using Sine Wave Electrostimulation at Varying | 178 files, 0.12 GB, species ['Rattus norvegicus - Norway rat'] 001033 0.240527.1432 | Characterization of Neuron Action Potential Triggering Using Sine Wave Electrostimulation at Varying | 193 files, 0.13 GB, species ['Rattus norvegicus - Norway rat'] 001034 0.240527.1431 | Characterization of Neuron Action Potential Triggering Using Sine Wave Electrostimulation at Varying | 178 files, 0.12 GB, species ['Rattus norvegicus - Norway rat'] 001035 0.240527.1432 | Characterization of Neuron Action Potential Triggering Using Sine Wave Electrostimulation at Varying | 254 files, 0.18 GB, species ['Rattus norvegicus - Norway rat']
001271 0.250805.1352 | Tetrodotoxin Blockade of Sodium Channels Abolishes Action Potential Formation During Intense Sinewav | 287 files, 0.2 GB, species ['Rattus norvegicus - Norway rat']
001272 0.250805.1358 | The Use of Di8-ANEPPS to Evaluate Neuronal Sensitivity to Unmodulated Sine Waves | 65 files, 0.04 GB, species ['Rattus norvegicus - Norway rat']
001273 0.250805.1356 | Kinetics of Membrane Charging and Its Dependence on Neuron Shape and Alignment with the Electric Fie | 181 files, 10.63 GB, species ['Rattus norvegicus - Norway rat']
001274 0.250805.1401 | Characterization of Plasma Membrane Polarization in Hippocampal Neurons During Sine Wave and Modulat | 232 files, 25.02 GB, species ['Rattus norvegicus - Norway rat']
001030: 17 files, 17 subject ids, species ['Mus musculus']
subject descriptions: ['Dissociated primary neurons from E18 C57BL/6 mouse hippocampal cells']
frames per file: {7000: 15, 6000: 2} | rates (frames/s): [3134.8]
stimuli: {('sine', 2000.0, nan): 17}
TTX sessions: 0 | labelled volts from 1.0 to 2.3
solution NaCl: 0 NaGluconate: 17
001032: 178 files, 14 subject ids, species ['Rattus norvegicus']
subject descriptions: ['Dissociated primary neurons from E18 Sprague Dawley rat hippocampal cells']
frames per file: {1200: 117, 900: 59, 1400: 1, 1800: 1} | rates (frames/s): [3086.4]
stimuli: {('mod', 2000.0, 20.0): 40, ('sine', 20.0, nan): 25, ('sine', 2000.0, nan): 80, ('sine', 20000.0, nan): 33}
TTX sessions: 0 | labelled volts from 2.5 to 18.0
solution NaCl: 178 NaGluconate: 0
001033: 193 files, 10 subject ids, species ['Rattus norvegicus']
subject descriptions: ['Dissociated primary neurons from E18 Sprague Dawley rat hippocampal cells']
frames per file: {1200: 193} | rates (frames/s): [3086.4]
stimuli: {('mod', 2000.0, 20.0): 31, ('mod', 20000.0, 20.0): 13, ('sine', 20.0, nan): 33, ('sine', 50.0, nan): 4, ('sine', 100.0, nan): 1, ('sine', 160.0, nan): 1, ('sine', 320.0, nan): 1, ('sine', 500.0, nan): 1, ('sine', 640.0, nan): 1, ('sine', 1280.0, nan): 1, ('sine', 2000.0, nan): 56, ('sine', 20000.0, nan): 47, ('sine', 40000.0, nan): 3}
TTX sessions: 0 | labelled volts from 2.0 to 20.0
solution NaCl: 193 NaGluconate: 0
001034: 178 files, 14 subject ids, species ['Rattus norvegicus']
subject descriptions: ['Dissociated primary neurons from E18 Sprague Dawley rat hippocampal cells']
frames per file: {1200: 117, 900: 59, 1800: 1, 1400: 1} | rates (frames/s): [3086.4]
stimuli: {('mod', 2000.0, 20.0): 40, ('sine', 20.0, nan): 25, ('sine', 2000.0, nan): 80, ('sine', 20000.0, nan): 33}
TTX sessions: 0 | labelled volts from 2.5 to 18.0
solution NaCl: 178 NaGluconate: 0
001035: 254 files, 15 subject ids, species ['Rattus norvegicus']
subject descriptions: ['Dissociated primary neurons from E18 C57BL/6 rat hippocampal cells']
frames per file: {1200: 253, 6000: 1} | rates (frames/s): [3086.4]
stimuli: {('sham', nan, nan): 1, ('sine', 20.0, nan): 63, ('sine', 40.0, nan): 3, ('sine', 80.0, nan): 4, ('sine', 160.0, nan): 4, ('sine', 320.0, nan): 2, ('sine', 640.0, nan): 3, ('sine', 1280.0, nan): 2, ('sine', 2000.0, nan): 73, ('sine', 2560.0, nan): 6, ('sine', 5120.0, nan): 5, ('sine', 10000.0, nan): 5, ('sine', 20000.0, nan): 81, ('sine', 100000.0, nan): 2}
TTX sessions: 0 | labelled volts from 2.0 to 20.0
solution NaCl: 254 NaGluconate: 0
001271: 287 files, 20 subject ids, species ['Rattus norvegicus']
subject descriptions: ['Dissociated primary neurons from E18 Sprague Dawley rat hippocampal cells']
frames per file: {1200: 287} | rates (frames/s): [3086.4]
stimuli: {('sham', nan, nan): 5, ('sine', 2000.0, nan): 282}
TTX sessions: 71 | labelled volts from 2.0 to 12.0
solution NaCl: 287 NaGluconate: 0
001272: 65 files, 3 subject ids, species ['Rattus norvegicus']
subject descriptions: ['Dissociated primary neurons from E18 Sprague Dawley rat hippocampal cells']
frames per file: {1200: 65} | rates (frames/s): [3086.4]
stimuli: {('sine', 2000.0, nan): 65}
TTX sessions: 0 | labelled volts from 2.0 to 14.0
solution NaCl: 65 NaGluconate: 0
# Recordings deposited more than once: identical pixel data (md5 of the whole stack).
h = pd.DataFrame([dict(dandiset=d, path=r["path"], md5=r["md5"]) for d in TL for r in REC[d] if "error" not in r])
dup = h[h.duplicated("md5", keep=False)].sort_values("md5")
print("files whose pixel data equal another file's:", len(dup), "in", dup.md5.nunique(), "groups")
print(dup.groupby("md5").dandiset.apply(lambda s: "+".join(sorted(s))).value_counts().to_dict())
pairs = dup.groupby("md5").path.apply(lambda s: tuple(sorted(x.split("/")[0] for x in s))).value_counts()
print(pairs.to_string())
print("Every 001034 file equals a 001032 file:", set(h[h.dandiset == "001034"].md5) <= set(h[h.dandiset == "001032"].md5))
files whose pixel data equal another file's: 356 in 178 groups
{'001032+001034': 178}
path
(sub-2-1, sub-2-3) 19
(sub-11-1, sub-11-3) 19
(sub-3-1, sub-3-3) 17
(sub-10-1, sub-10-3) 17
(sub-14-1, sub-14-3) 15
(sub-5-1, sub-5-3) 15
(sub-7-1, sub-7-3) 15
(sub-6-1, sub-6-3) 14
(sub-13-1, sub-13-3) 13
(sub-9-1, sub-9-3) 10
(sub-1-1, sub-1-3) 9
(sub-12-1, sub-12-3) 7
(sub-4-1, sub-4-3) 4
(sub-8-1, sub-8-3) 4
Every 001034 file equals a 001032 file: True
# Units and ranges of the measured variable (soma dF, %), the ROI sizes, and the pre-stimulus noise.
T = pd.concat([TAB[d] for d in TL if d != "001034"], ignore_index=True)
T1 = T[(T.o3 == "a")]
print(T1.groupby(["dandiset", "o1"]).agg(sessions=("path", "nunique"), roi_pixels_median=("npix", "median"), roi_pixels_min=("npix", "min"),
pre_sd_median=("pre_sd", "median")).round(2).to_string())
print("dF range over all sessions, ROI c, F0 a:")
for d in TL:
if d == "001034": continue
v = np.concatenate([r["traces"][("c", "a")][0] for r in REC[d] if "error" not in r])
print(f" {d}: {v.min():.1f} to {v.max():.1f} % (1st/99th percentile {np.percentile(v, 1):.1f} / {np.percentile(v, 99):.1f})")
f0 = [r["frame0_ratio"] for d in TL if d != "001034" for r in REC[d] if "error" not in r]
print("frame 0 over frame 1 (mean of the image), every file: min %.2f, median %.2f; frame 0 is dropped" % (min(f0), np.median(f0)))
sessions roi_pixels_median roi_pixels_min pre_sd_median
dandiset o1
001030 a 17 15.0 1 2.49
b 17 8.0 1 3.43
c 17 256.0 256 1.11
001032 a 178 36.0 2 1.71
b 178 10.0 2 3.10
c 178 256.0 256 0.93
001033 a 193 39.0 1 1.58
b 193 14.0 1 2.80
c 193 256.0 256 0.81
001035 a 254 35.0 1 1.57
b 254 6.0 1 3.50
c 254 256.0 256 0.79
001271 a 287 34.0 1 1.54
b 287 13.0 1 2.51
c 287 256.0 256 0.81
001272 a 65 31.0 26 1.18
b 65 12.0 4 1.56
c 65 256.0 256 0.55
dF range over all sessions, ROI c, F0 a:
001030: -8.6 to 16.9 % (1st/99th percentile -6.3 / 9.3)
001032: -89.7 to 24.3 % (1st/99th percentile -2.7 / 11.5)
001033: -5.7 to 16.4 % (1st/99th percentile -2.7 / 8.5)
001035: -12.8 to 16.1 % (1st/99th percentile -3.9 / 8.1)
001271: -5.8 to 17.5 % (1st/99th percentile -2.5 / 8.2)
001272: -5.3 to 5.9 % (1st/99th percentile -0.9 / 3.9)
frame 0 over frame 1 (mean of the image), every file: min 1.51, median 3.15; frame 0 is dropped
# The clock: are APs locked to the stated stimulus onset? AP rate per 10-ms bin, before and after onset,
# over all 250-ms sessions with a stimulus (ROI c, F0 a, O3 a), and in sham sessions. APs are detected
# over the whole recording here (not only in the stimulus window) to show the contrast.
def all_peaks(d):
out = []
for r in REC[d]:
if "error" in r or r["stim"]["kind"] not in ("sine", "mod", "sham"): continue
dF, info = r["traces"][("c", "a")]
pt, ph, pw = ts.detect_aps(r["t"], dF, 0, 1e9, info["pre_sd"], "a")
out.append((r["stim"]["kind"], pt - r["onset"], r["t"][-1] - r["onset"], r["onset"]))
return out
for d in ["001032", "001033", "001035", "001271", "001030"]:
P = all_peaks(d)
stim = [p for p in P if p[0] != "sham"]
pre = sum(((x[1] > -x[3] + 3) & (x[1] < 0)).sum() for x in stim); pre_t = sum(x[3] - 3 for x in stim) / 1000
dur = 2000 if d == "001030" else 250
dur_n = sum(((x[1] >= 0) & (x[1] < dur)).sum() for x in stim); dur_t = sum(min(dur, x[2]) for x in stim) / 1000
post = sum(((x[1] >= dur + 20)).sum() for x in stim); post_t = sum(max(0, x[2] - dur - 20) for x in stim) / 1000
print(f"{d}: AP rate before onset {pre / pre_t:.2f}/s, during stimulus {dur_n / dur_t:.2f}/s, from 20 ms after its end {post / max(post_t, 1e-9):.2f}/s "
f"(sessions {len(stim)}); sham sessions {sum(p[0] == 'sham' for p in P)}, APs in them {sum(len(p[1]) for p in P if p[0] == 'sham')}")
if d == "001033":
lat = np.concatenate([x[1] for x in stim]); h_, e_ = np.histogram(lat, bins=np.arange(-50, 341, 10))
print(" 001033 AP counts per 10 ms from -50 ms:", h_.tolist())
001032: AP rate before onset 0.00/s, during stimulus 15.00/s, from 20 ms after its end 0.36/s (sessions 178); sham sessions 0, APs in them 0 001033: AP rate before onset 0.00/s, during stimulus 12.62/s, from 20 ms after its end 0.45/s (sessions 193); sham sessions 0, APs in them 0 001033 AP counts per 10 ms from -50 ms: [0, 0, 0, 0, 0, 42, 43, 28, 28, 27, 25, 15, 27, 18, 20, 16, 20, 34, 31, 10, 19, 16, 33, 22, 22, 21, 22, 23, 28, 19, 5, 0, 1, 1, 1, 1, 0, 1, 1] 001035: AP rate before onset 0.00/s, during stimulus 12.35/s, from 20 ms after its end 0.06/s (sessions 253); sham sessions 1, APs in them 0 001271: AP rate before onset 0.00/s, during stimulus 18.62/s, from 20 ms after its end 0.00/s (sessions 282); sham sessions 5, APs in them 0 001030: AP rate before onset 0.00/s, during stimulus 16.28/s, from 20 ms after its end 0.00/s (sessions 17); sham sessions 0, APs in them 0
import re
# The TTX and vehicle wording of the session descriptions (round-2 verifier: the findings quote it).
for d in ["001271", "001272"]:
R = [r for r in REC[d] if "error" not in r]
for pat in ["caused sustained depolarization", "prevented", "vehicle control"]:
print(d, repr(pat), "in", sum(pat in r["session_description"] for r in R), "of", len(R), "session descriptions")
for e in sorted({m.strip() for r in R for m in re.findall(r"[^.]*(?:sustained depolarization|vehicle control)[^.]*\.", r["session_description"])})[:4]:
print(" ", e[:300])
for d in ["001271", "001272"]:
desc = ts.dandiset_meta(d).get("description", "")
print(d, "dandiset description, sentences on depolarization or vehicle:",
[e.strip()[:300] for e in re.findall(r"[^.]*(?:depolarization|vehicle)[^.]*\.", desc)])
001271 'caused sustained depolarization' in 71 of 287 session descriptions
001271 'prevented' in 0 of 287 session descriptions
001271 'vehicle control' in 216 of 287 session descriptions
The experiment was performed in the presence of the sodium channel inhibitor 1 µM tetrodotoxin, which abolished the action potentials and caused sustained depolarization during application of stimulating waves.
The experiment was performed without tetrodotoxin (vehicle control).
001272 'caused sustained depolarization' in 0 of 65 session descriptions
001272 'prevented' in 0 of 65 session descriptions
001272 'vehicle control' in 65 of 65 session descriptions
The experiment was performed without tetrodotoxin (vehicle control).
001271 dandiset description, sentences on depolarization or vehicle: ['Applying 1 µM TTX suppressed action potentials, as expected, and completely prevented sustained depolarization.']
001272 dandiset description, sentences on depolarization or vehicle: []
Checks of the metadata against the data:
- 001048 holds TIRF images of HEK cells under patch clamp, with species human; it is another study.
Its only link to this paper is the related-resource entry. The paper's data are in the ten
dandisets its Data and code availability lists, and the claims file marks every claim
data: elsewherein them. - 001032 and 001034 hold the same 178 recordings, pixel for pixel, under different subject names (sub-N-1 in 001032 and sub-N-3 in 001034). They are counted once, from 001032.
- The images are 16 x 16 pixels, where the paper describes 64 x 64 frames after binning; the depositors cropped each recording around the soma. Frame 0 of every file is a bright artefact and is dropped.
- The frame rates are 3,086.4 frames/s (001032 to 001035, 001271, 001272) and 3,134.8 (001030), close to the paper's 3,086.4 and 3,130.
- The APs start after the stated onset (50 ms; 25 ms in 001030) and the rate before it is low: this is consistent with the stated stimulus timing, which is not recorded in the files (no stimulus series). Sham sessions have no APs.
- 001030 states mouse (Mus musculus, C57BL/6) neurons; the paper used only rat neurons. 001035 calls its neurons "C57BL/6 rat". 001030 makes each session a subject (sub-c7-1-45V); the neuron is read from the name.
- The voltages are labelled "peak-to-peak" in 001032 to 001035 and "one polarity" in 001030, and the files say they are "not the actual amplitude measured on the oscilloscope", with frequency dependent calibration "reported separately". The calibration is not in the dandisets.
- There are no NaGlu (low-conductivity) sessions with 250-ms waves anywhere: Figures 3A and 3B are not deposited. 001030 holds only the 2-s ramped waves of Figure 3C.
- 001272 holds 3 Di-8-ANEPPS neurons; the paper reports 8.
- 001032, 001033 and 001035 also hold sine waves of 20 Hz to 10 kHz, and 001033 and 001035 a few at 40 and 100 kHz, that the paper does not describe ("20Hz sine" in most neurons, frequency series in a few). They are not used.
- The TTX session descriptions of 001271 (71 sessions) say TTX "abolished the action potentials and caused sustained depolarization", where the dandiset's description says TTX "completely prevented sustained depolarization". Every session description of 001272 says "without tetrodotoxin (vehicle control)", the wording of the 001271 vehicle sessions.
Readings and selections used by the claims¶
All claims on the time-lapse data are run under every combination of the readings the claims file
names: ROI (O1) (a) the Otsu component nearest the centre, (b) that mask eroded by 2 pixels; F0 (O2)
(a) all frames before onset, (b) the last 20 ms; AP height (O3) (a) above 5% dF, (b) above 5 SD of
the pre-onset dF. A third ROI reading (c) is added: the whole 16 x 16 crop that the depositors cut
around the soma. The paper measures "the emission intensity over the entire cell soma" (Methods,
Recording of the optical APs); the crop is the dataset's own region of the soma, and the claims file
says ROIs are taken from the dataset when it holds them (O1). The Otsu masks are small (a median of
about 35 of 256 pixels, some of 1 to 3 pixels, above) and their pre-onset noise is about twice the
crop's. The claims file's header blocks a claim whose readings disagree; protocol step 1.6 allows
adopting the best-supported choice as an assumption instead. Where the ROI readings give different
outcomes, the outcome is the one under the crop ROI, recorded as resting on that assumption
(assumed), and every reading is shown.
AP detection (O3): a local maximum of dF above the height reading, at least 5 ms from a larger one, with a width at half height of 2 to 20 ms. The width is measured both at half the peak's dF (its height above F0) and at half its prominence (its rise above the neighbouring troughs), and must be 2 to 20 ms under both. Each measure alone fails here: at half the dF, small wiggles on the tail of an AP count as APs 5 ms later (maximum rates of 136 to 182 Hz and up to 18 APs in a modulated session); at half the prominence, small oscillations on a raised plateau count as APs, so that the claims file's sustained depolarization (no AP in the last 100 ms) can hardly occur. Requiring both is an assumption made here, in line with the paper's "brief emission upticks"; both single measures are run at the end (AP detector variants).
The 250-ms FluoVolt neurons are those of 001032, 001033 and 001035 (001034 is a copy of 001032). Stimuli are 2 and 20 kHz, unmodulated ("sine") or modulated at 20 Hz ("spindle" in the files). Traces are not corrected for bleaching: the claims file (F2) says the paper corrected only the 2-s recordings of Figure 3C.
RD = [(o1, o2, o3) for o1 in "cab" for o2 in "ab" for o3 in "ab"]
FL_ALL = pd.concat([TAB[d] for d in ["001032", "001033", "001035"]], ignore_index=True)
FL = FL_ALL[FL_ALL.kind.isin(["sine", "mod"]) & FL_ALL.carrier.isin([2000.0, 20000.0])].copy()
FL["stim"] = np.where(FL["kind"] == "mod", "mod", "pure") + "_" + (FL.carrier / 1000).astype(int).astype(str) + "k"
def rd(df, r): return df[(df.o1 == r[0]) & (df.o2 == r[1]) & (df.o3 == r[2])]
g0 = rd(FL, ("c", "a", "a"))
print("neurons with 2/20 kHz sessions:", g0.neuron.nunique(), "| sessions:", len(g0))
print(g0.groupby("stim").agg(neurons=("neuron", "nunique"), sessions=("path", "size")).to_string())
print("by dandiset:", g0.groupby("dandiset").neuron.nunique().to_dict())
print("sessions ending before the stimulus ends (900 frames, 291 ms):", int((g0.t_end < 300).sum()))
neurons with 2/20 kHz sessions: 34 | sessions: 454
neurons sessions
stim
mod_20k 3 13
mod_2k 14 71
pure_20k 17 161
pure_2k 28 209
by dandiset: {'001032': 13, '001033': 7, '001035': 14}
sessions ending before the stimulus ends (900 frames, 291 ms): 38
ap-selection¶
Neurons were selected for APs peaking above 5% dF and lasting about 10 ms; typical APs were 11.3% dF high and 8 ms long, excluding the afterhyperpolarization (mean of 22).
Per neuron, the mean peak dF and duration of its APs in the 250-ms FluoVolt sessions (every stimulus in 001032, 001033 and 001035), then the mean and SEM over neurons. Which APs: (a) every AP, (b) the first AP of each session. Duration: (a) width at half height, (b) from the rise to the first return to F0. Expected: 11.3 +- 0.9 % dF and 8 +- 0.7 ms.
rows = []
for r in RD:
g = rd(FL_ALL[FL_ALL.kind.isin(["sine", "mod"])], r)
g = g[g.n_ap > 0]
for which in ["every", "first"]:
for dur in ["half_width", "return"]:
per = []
for n, x in g.groupby("neuron"):
if which == "every":
h = np.concatenate(x.ap_h.values); w = np.concatenate((x.ap_w if dur == "half_width" else x.ap_ret).values)
else:
h = np.array([v[0] for v in x.ap_h]); w = np.array([v[0] for v in (x.ap_w if dur == "half_width" else x.ap_ret)])
per.append((h.mean(), w.mean()))
per = np.array(per)
rows.append(dict(o1=r[0], o2=r[1], o3=r[2], aps=which, duration=dur, neurons=len(per),
peak=per[:, 0].mean(), peak_sem=stats.sem(per[:, 0]), dur_ms=per[:, 1].mean(), dur_sem=stats.sem(per[:, 1])))
APS = pd.DataFrame(rows)
APS["peak_ok"] = (APS.peak - 11.3).abs() <= 0.9
APS["dur_ok"] = (APS.dur_ms - 8).abs() <= 0.7
print(APS.round(2).to_string())
print("\nreadings with both within the SEM:", int((APS.peak_ok & APS.dur_ok).sum()), "of", len(APS))
for o1 in "cab":
x = APS[APS.o1 == o1]
print(f"ROI {o1}: peak within the SEM in {int(x.peak_ok.sum())} of {len(x)} readings, duration in {int(x.dur_ok.sum())}; "
f"half width {x[x.duration == 'half_width'].dur_ms.min():.2f} to {x[x.duration == 'half_width'].dur_ms.max():.2f} ms, "
f"rise to return {x[x.duration == 'return'].dur_ms.min():.2f} to {x[x.duration == 'return'].dur_ms.max():.2f} ms")
print("ROI c: peak %.2f to %.2f %%; duration %.2f to %.2f ms" % (APS[APS.o1 == "c"].peak.min(), APS[APS.o1 == "c"].peak.max(),
APS[APS.o1 == "c"].dur_ms.min(), APS[APS.o1 == "c"].dur_ms.max()))
o1 o2 o3 aps duration neurons peak peak_sem dur_ms dur_sem peak_ok dur_ok 0 c a a every half_width 39 10.53 0.38 3.70 0.18 True False 1 c a a every return 39 10.53 0.38 23.51 4.15 True False 2 c a a first half_width 39 11.46 0.41 3.16 0.11 True False 3 c a a first return 39 11.46 0.41 18.12 2.42 True False 4 c a b every half_width 39 10.39 0.39 3.71 0.18 False False 5 c a b every return 39 10.39 0.39 24.14 4.26 False False 6 c a b first half_width 39 11.29 0.43 3.15 0.11 True False 7 c a b first return 39 11.29 0.43 18.05 2.42 True False 8 c b a every half_width 39 10.41 0.37 3.65 0.18 True False 9 c b a every return 39 10.41 0.37 23.05 4.09 True False 10 c b a first half_width 39 11.32 0.41 3.12 0.11 True False 11 c b a first return 39 11.32 0.41 17.38 2.30 True False 12 c b b every half_width 39 10.27 0.39 3.66 0.17 False False 13 c b b every return 39 10.27 0.39 23.65 4.21 False False 14 c b b first half_width 39 11.16 0.43 3.11 0.11 True False 15 c b b first return 39 11.16 0.43 17.31 2.31 True False 16 a a a every half_width 37 11.18 0.71 2.97 0.12 True False 17 a a a every return 37 11.18 0.71 16.13 4.21 True False 18 a a a first half_width 37 11.55 0.71 2.81 0.09 True False 19 a a a first return 37 11.55 0.71 12.11 2.65 True False 20 a a b every half_width 29 11.21 0.61 3.05 0.14 True False 21 a a b every return 29 11.21 0.61 18.68 5.30 True False 22 a a b first half_width 29 11.56 0.64 2.85 0.11 True False 23 a a b first return 29 11.56 0.64 13.71 3.34 True False 24 a b a every half_width 37 11.07 0.70 2.92 0.11 True False 25 a b a every return 37 11.07 0.70 15.54 4.09 True False 26 a b a first half_width 37 11.41 0.70 2.76 0.08 True False 27 a b a first return 37 11.41 0.70 11.95 2.61 True False 28 a b b every half_width 29 11.01 0.58 3.00 0.13 True False 29 a b b every return 29 11.01 0.58 18.06 5.15 True False 30 a b b first half_width 29 11.33 0.62 2.81 0.10 True False 31 a b b first return 29 11.33 0.62 13.70 3.29 True False 32 b a a every half_width 34 12.79 0.67 2.60 0.06 False False 33 b a a every return 34 12.79 0.67 10.77 2.95 False False 34 b a a first half_width 34 12.82 0.67 2.54 0.06 False False 35 b a a first return 34 12.82 0.67 9.20 2.21 False False 36 b a b every half_width 18 14.42 0.65 2.75 0.12 False False 37 b a b every return 18 14.42 0.65 10.51 1.44 False False 38 b a b first half_width 18 14.65 0.67 2.64 0.12 False False 39 b a b first return 18 14.65 0.67 8.57 0.87 False True 40 b b a every half_width 34 12.68 0.66 2.59 0.06 False False 41 b b a every return 34 12.68 0.66 8.18 0.88 False True 42 b b a first half_width 34 12.69 0.66 2.52 0.06 False False 43 b b a first return 34 12.69 0.66 7.30 0.67 False True 44 b b b every half_width 19 14.41 0.61 2.71 0.12 False False 45 b b b every return 19 14.41 0.61 10.78 1.39 False False 46 b b b first half_width 19 14.67 0.62 2.60 0.12 False False 47 b b b first return 19 14.67 0.62 9.08 0.93 False False readings with both within the SEM: 0 of 48 ROI c: peak within the SEM in 12 of 16 readings, duration in 0; half width 3.11 to 3.71 ms, rise to return 17.31 to 24.14 ms ROI a: peak within the SEM in 16 of 16 readings, duration in 0; half width 2.76 to 3.05 ms, rise to return 11.95 to 18.68 ms ROI b: peak within the SEM in 0 of 16 readings, duration in 3; half width 2.52 to 2.75 ms, rise to return 7.30 to 10.78 ms ROI c: peak 10.27 to 11.46 %; duration 3.11 to 24.14 ms
# Per dandiset (the neurons come from three dandisets), ROI c, F0 a, O3 a, every AP.
g = rd(FL_ALL[FL_ALL.kind.isin(["sine", "mod"])], ("c", "a", "a")); g = g[g.n_ap > 0]
for d, x in g.groupby("dandiset"):
per = x.groupby("neuron").apply(lambda y: pd.Series(dict(h=np.concatenate(y.ap_h.values).mean(), w=np.concatenate(y.ap_w.values).mean(),
ret=np.concatenate(y.ap_ret.values).mean())), include_groups=False)
print(d, len(per), "neurons: peak %.1f, half width %.1f ms, rise to return %.1f ms" % tuple(per.mean()))
001032 14 neurons: peak 12.1, half width 4.3 ms, rise to return 30.6 ms 001033 10 neurons: peak 9.3, half width 3.5 ms, rise to return 23.5 ms 001035 15 neurons: peak 9.9, half width 3.3 ms, rise to return 16.9 ms
Result: partially reproduced. Over 39 neurons the mean AP peak is within 11.3 +- 0.9 % dF in 12 of 16 readings under the crop ROI (c), 10.27 to 11.46 %, but the duration matches 8 +- 0.7 ms under no reading with that ROI: the width at half height is 3.11 to 3.71 ms and the rise-to-return duration 17.31 to 24.14 ms, the latter lengthened by the noise around F0. Under ROI (a) the peak is within the SEM in all 16 readings (11.01 to 11.56 %) and the duration in none (2.76 to 3.05 ms half width, 11.97 to 18.68 ms rise to return). Under ROI (b) the peak is 12.68 to 14.67 %, outside it, and three readings give a rise-to-return duration in range (7.32 to 8.57 ms). The peak holds under the crop ROI, taken as the paper's whole soma (an assumption); the duration holds under no reading in which the peak holds.
fig1-pure-sine-excites-all¶
Unmodulated 250-ms sine waves at 2 and 20 kHz excited every neuron tested; near threshold isolated APs came 100 to 200 ms after onset, and stronger fields shortened the delay and evoked sustained firing whose maximum rate rose linearly to a ceiling of 60 to 90 Hz.
Per neuron and frequency (unmodulated sessions): excited when any session has an AP; the threshold is the lowest labelled voltage with an AP (F5) and the latency there the mean first-AP latency of the sessions at that voltage. Spearman rho of latency on voltage over a neuron's sessions with APs (at least 3 sessions at 2 or more voltages, a standard choice). The highest maximum firing rate of each neuron (inverse of its first inter-AP interval, F4) over both frequencies, readings (a) the median over neurons in 54 to 99 Hz, (b) at least 75% of neurons in 54 to 99 Hz.
def fig1(g):
g = g[g.kind == "sine"]
out = {}
for f, x in g.groupby("carrier"):
ex = x.groupby("neuron").n_ap.max() > 0
out[f"excited_{int(f/1000)}k"] = ex.mean() * 100; out[f"n_{int(f/1000)}k"] = len(ex)
lat, rho = [], []
for (n, f), x in g.groupby(["neuron", "carrier"]):
y = x[x.n_ap > 0]
if len(y) == 0: continue
v0 = y.volts.min(); lat.append(y[y.volts == v0].latency.mean())
if len(y) >= 3 and y.volts.nunique() >= 2:
rho.append(stats.spearmanr(y.volts, y.latency)[0])
rho = np.array(rho); lat = np.array(lat)
mr = g.groupby("neuron").max_rate.max().dropna()
out.update(median_latency_at_threshold=np.median(lat), n_lat=len(lat), rho_neg=(rho < 0).sum(), n_rho=len(rho),
rho_median=np.median(rho), max_rate_median=mr.median(), frac_rate_in_54_99=((mr >= 54) & (mr <= 99)).mean() * 100,
n_rate=len(mr), max_rate_min=mr.min(), max_rate_max=mr.max())
out["holds_a"] = (out["excited_2k"] >= 90 and out["excited_20k"] >= 90 and 90 <= out["median_latency_at_threshold"] <= 220
and out["rho_neg"] > out["n_rho"] / 2 and 54 <= out["max_rate_median"] <= 99)
out["holds_b"] = (out["excited_2k"] >= 90 and out["excited_20k"] >= 90 and 90 <= out["median_latency_at_threshold"] <= 220
and out["rho_neg"] > out["n_rho"] / 2 and out["frac_rate_in_54_99"] >= 75)
return out
F1 = pd.DataFrame([dict(o1=r[0], o2=r[1], o3=r[2], **fig1(rd(FL, r))) for r in RD])
print(F1.round(1).to_string())
o1 o2 o3 excited_2k n_2k excited_20k n_20k median_latency_at_threshold n_lat rho_neg n_rho rho_median max_rate_median frac_rate_in_54_99 n_rate max_rate_min max_rate_max holds_a holds_b 0 c a a 96.4 28 88.2 17 60.8 42 34 37 -0.9 41.7 31.0 29 5.9 154.3 False False 1 c a b 96.4 28 94.1 17 59.8 43 35 38 -0.9 40.1 30.0 30 5.9 154.3 False False 2 c b a 96.4 28 88.2 17 63.9 42 34 37 -0.9 42.3 28.6 28 16.1 154.3 False False 3 c b b 96.4 28 94.1 17 61.8 43 35 38 -0.9 41.7 27.6 29 16.1 154.3 False False 4 a a a 85.7 28 58.8 17 64.4 34 23 29 -0.7 32.0 9.1 22 4.3 154.3 False False 5 a a b 71.4 28 58.8 17 64.4 30 21 24 -0.8 32.8 10.5 19 4.3 154.3 False False 6 a b a 85.7 28 58.8 17 88.4 34 23 29 -0.7 30.6 9.1 22 4.3 154.3 False False 7 a b b 67.9 28 58.8 17 66.6 29 20 24 -0.8 32.8 10.5 19 4.3 154.3 False False 8 b a a 78.6 28 52.9 17 61.5 31 10 16 -0.4 23.0 11.8 17 10.2 85.7 False False 9 b a b 42.9 28 23.5 17 76.2 16 6 7 -0.5 22.6 20.0 10 9.1 85.7 False False 10 b b a 78.6 28 52.9 17 61.5 31 10 16 -0.4 23.0 11.8 17 10.2 85.7 False False 11 b b b 39.3 28 23.5 17 78.6 15 5 6 -0.6 22.6 20.0 10 9.1 85.7 False False
# Per frequency (ROI c, F0 a, O3 a), the per-dandiset breakdown, and the linear rise of the maximum
# rate below each neuron's ceiling (reported, does not decide).
g = rd(FL, ("c", "a", "a")); g = g[g.kind == "sine"]
for f, x in g.groupby("carrier"):
lat = [y[y.volts == y.volts.min()].latency.mean() for n, y in x[x.n_ap > 0].groupby("neuron")]
mr = x.groupby("neuron").max_rate.max().dropna()
print(f"{int(f/1000)} kHz: {x.neuron.nunique()} neurons, excited {x.groupby('neuron').n_ap.max().gt(0).sum()}; median latency at threshold "
f"{np.median(lat):.0f} ms (IQR {np.percentile(lat, 25):.0f} to {np.percentile(lat, 75):.0f}); median highest max rate {mr.median():.0f} Hz")
for d, x in g.groupby("dandiset"):
print(d, "neurons", x.neuron.nunique(), "excited", x.groupby("neuron").n_ap.max().gt(0).sum(),
"median highest max rate %.0f Hz" % x.groupby("neuron").max_rate.max().median())
sl = []
for n, x in g.groupby("neuron"):
y = x.dropna(subset=["max_rate"])
if len(y) < 3: continue
vmax = y.loc[y.max_rate.idxmax(), "volts"]; y = y[y.volts <= vmax]
if y.volts.nunique() >= 3: sl.append(stats.linregress(y.volts, y.max_rate).slope)
sl = np.array(sl)
print(f"linear fit of max rate on voltage up to the ceiling: {len(sl)} neurons, slope median {np.median(sl):.1f} Hz/V, positive in {(sl > 0).sum()}")
fig, ax = plt.subplots(1, 2, figsize=(10, 3.5))
for n, x in g.groupby("neuron"):
y = x.dropna(subset=["max_rate"]).sort_values("volts")
ax[0].plot(y.volts, y.max_rate, ".-", lw=.5, ms=3)
y = x[x.n_ap > 0].sort_values("volts"); ax[1].plot(y.volts, y.latency, ".-", lw=.5, ms=3)
ax[0].set(xlabel="labelled volts (peak to peak)", ylabel="max firing rate (Hz)"); ax[1].set(xlabel="labelled volts", ylabel="first-AP latency (ms)")
plt.tight_layout(); plt.show()
print("plotted max rate range %.0f to %.0f Hz; latency range %.0f to %.0f ms; nothing outside the axes"
% (g.max_rate.min(), g.max_rate.max(), g.latency.min(), g.latency.max()))
2 kHz: 28 neurons, excited 27; median latency at threshold 59 ms (IQR 34 to 115); median highest max rate 31 Hz 20 kHz: 17 neurons, excited 15; median latency at threshold 124 ms (IQR 48 to 165); median highest max rate 32 Hz 001032 neurons 11 excited 11 median highest max rate 37 Hz 001033 neurons 7 excited 7 median highest max rate 39 Hz 001035 neurons 14 excited 13 median highest max rate 45 Hz linear fit of max rate on voltage up to the ceiling: 23 neurons, slope median 12.3 Hz/V, positive in 23
plotted max rate range 5 to 154 Hz; latency range 1 to 242 ms; nothing outside the axes
# First APs within 5 ms of the stimulus onset (round-1 verifier, problem 7): are they artefacts of the
# onset? The same stimulus under TTX, which removes the APs, shows no such transient: the maximum dF in
# the 10 ms after onset, in the 001271 sessions with and without TTX (ROI c, F0 a).
g = rd(FL, ("c", "a", "a")); g = g[g.kind.isin(["sine", "mod"])]
e = g[g.latency < 5]
print(f"sessions with a first AP within 5 ms of onset: {len(e)} of {int((g.n_ap > 0).sum())} with APs, in {e.neuron.nunique()} neurons; "
f"labelled volts median {e.volts.median():.1f} (sessions with APs: {g[g.n_ap > 0].volts.median():.1f}); APs per session median {e.n_ap.median():.0f} "
f"(all sessions with APs: {g[g.n_ap > 0].n_ap.median():.0f}); first-AP height median {np.median([h[0] for h in e.ap_h]):.1f} %, half width {np.median([w[0] for w in e.ap_w]):.1f} ms")
for ttx_ in [True, False]:
mx = []
for r in REC["001271"]:
if "error" in r or r["stim"]["kind"] != "sine" or r["stim"]["ttx"] != ttx_: continue
dF, _ = r["traces"][("c", "a")]; t = r["t"]
mx.append(float(dF[(t >= 50) & (t < 60)].max()))
mx = np.array(mx)
print(f"001271 {'with' if ttx_ else 'without'} TTX: {len(mx)} sessions, max dF 0 to 10 ms after onset median {np.median(mx):.2f} %, 95th percentile {np.percentile(mx, 95):.2f} %, above 5 %: {int((mx > 5).sum())}")
th1 = []
for (n, f), x in g[g.kind == "sine"].groupby(["neuron", "carrier"]):
y = x[x.n_ap > 0]
if len(y): th1.append(y[y.volts == y.volts.min()].latency.mean())
th1 = np.array(th1)
print(f"latency at threshold: {len(th1)} neuron-frequencies, {int((th1 < 5).sum())} under 5 ms; median {np.median(th1):.1f} ms, without those {np.median(th1[th1 >= 5]):.1f} ms")
sessions with a first AP within 5 ms of onset: 66 of 325 with APs, in 16 neurons; labelled volts median 7.0 (sessions with APs: 6.5); APs per session median 10 (all sessions with APs: 4); first-AP height median 12.4 %, half width 3.1 ms 001271 with TTX: 71 sessions, max dF 0 to 10 ms after onset median 1.90 %, 95th percentile 3.51 %, above 5 %: 0 001271 without TTX: 211 sessions, max dF 0 to 10 ms after onset median 7.93 %, 95th percentile 13.80 %, above 5 %: 145
latency at threshold: 42 neuron-frequencies, 2 under 5 ms; median 60.8 ms, without those 63.9 ms
Result: partially reproduced. Under the crop ROI unmodulated waves excited 27 of 28 neurons at 2 kHz (96.4 %) and 15 of 17 at 20 kHz (88.2 % with O3 a, 94.1 % with O3 b), and the first-AP latency fell with voltage in 34 of 37 neuron-frequencies (35 of 38 with O3 b). The other two parts fail under every reading: the median latency at threshold is 59.8 to 88.4 ms (59 ms at 2 kHz and 124 ms at 20 kHz under ROI c), and the median highest maximum firing rate is 22.6 to 42.3 Hz, with 9.1 to 31.0 % of neurons in 54 to 99 Hz. Two of the 42 latencies at threshold (ROI c) are under 5 ms; without them the median is 63.9 ms instead of 60.8. APs this early are not onset artefacts: the same stimulus under TTX gives no transient above 5 % in the first 10 ms in any of 71 sessions, against 145 of 211 without TTX. Under the Otsu ROIs fewer neurons reach the AP height (39.3 to 85.7 % excited at 2 kHz). Below each neuron's highest rate the rate rises with voltage (slope positive in 23 of 23 neurons, not deciding).
fig1-excitation-block¶
At the highest fields, unmodulated sine waves merged APs into low-amplitude oscillations on a continuous depolarization, which could flatten and leave the membrane depolarized for the rest of the stimulus (an excitation block), and this depolarization reversed immediately when the stimulus stopped.
Per neuron and frequency (unmodulated): the voltage giving the highest maximum firing rate; neurons tested above it; sustained depolarization (O4: no AP in the last 100 ms of the stimulus and mean dF there above (a) 5 SD of pre-onset dF or (b) half the neuron's median AP height) at their highest voltage; and the mean dF over the 50 ms after the stimulus, below the O4 level. Sessions of 900 frames end at 291 ms, before the stimulus ends, so their post-stimulus window is missing.
def block(g, o4):
g = g[g.kind == "sine"]
medh = g.groupby("neuron").ap_h.apply(lambda s: np.median(np.concatenate(s.values)) if len(np.concatenate(s.values)) else np.nan)
res = []
for (n, f), x in g.groupby(["neuron", "carrier"]):
y = x.dropna(subset=["max_rate"])
if len(y) == 0: continue
vbest = y.loc[y.max_rate.idxmax(), "volts"]
if x.volts.max() <= vbest: continue
top = x[x.volts == x.volts.max()]
s = [ts.sustained(row, medh[n], o4) for _, row in top.iterrows()]
sus = any(v[0] for v in s)
post = [row.mean_post50 < lvl for (_, row), (ok, lvl) in zip(top.iterrows(), s) if ok and np.isfinite(row.mean_post50)]
res.append(dict(neuron=n, carrier=f, sustained=sus, post_below=all(post) if post else np.nan, n_post=len(post)))
return pd.DataFrame(res)
rows = []
for r in RD:
for o4 in "ab":
b = block(rd(FL, r), o4)
rows.append(dict(o1=r[0], o2=r[1], o3=r[2], o4=o4, tested_above=len(b), sustained=int(b.sustained.sum()),
post_checked=int(b.n_post.gt(0).sum()), post_below=int((b.post_below == True).sum()),
holds=b.sustained.sum() > len(b) / 2 and (b.post_below == True).sum() == b.n_post.gt(0).sum()))
BLK = pd.DataFrame(rows)
print(BLK.to_string())
o1 o2 o3 o4 tested_above sustained post_checked post_below holds 0 c a a a 14 0 0 0 False 1 c a a b 14 0 0 0 False 2 c a b a 14 0 0 0 False 3 c a b b 14 0 0 0 False 4 c b a a 12 0 0 0 False 5 c b a b 12 0 0 0 False 6 c b b a 12 0 0 0 False 7 c b b b 12 0 0 0 False 8 a a a a 15 2 2 1 False 9 a a a b 15 2 2 1 False 10 a a b a 12 2 2 1 False 11 a a b b 12 2 2 1 False 12 a b a a 14 2 2 1 False 13 a b a b 14 2 2 1 False 14 a b b a 13 2 2 1 False 15 a b b b 13 2 2 1 False 16 b a a a 13 0 0 0 False 17 b a a b 13 1 1 1 False 18 b a b a 4 0 0 0 False 19 b a b b 4 1 1 1 False 20 b b a a 13 0 0 0 False 21 b b a b 13 1 1 1 False 22 b b b a 4 0 0 0 False 23 b b b b 4 1 1 1 False
# Mean dF in the last 100 ms of the stimulus at each neuron's highest voltage, and APs there (ROI c, F0 a, O3 a).
g = rd(FL, ("c", "a", "a")); g = g[g.kind == "sine"]
top = g.loc[g.groupby(["neuron", "carrier"]).volts.transform("max") == g.volts]
print(top[["neuron", "carrier", "volts", "n_ap", "n_ap_last100", "mean_last100", "pre_sd", "mean_post50"]].round(2).to_string())
# Variant (does not decide): sessions at any voltage that meet O4 (a), and the rise of the late mean dF with voltage.
g2 = g.assign(sus=(g.n_ap_last100 == 0) & (g.mean_last100 > 5 * g.pre_sd))
print("sessions meeting O4 (a) at any voltage:", int(g2.sus.sum()), "in", g2[g2.sus].neuron.nunique(), "of", g.neuron.nunique(), "neurons")
print(g2[g2.sus][["neuron", "carrier", "volts", "n_ap", "mean_last100", "pre_sd", "t_end"]].round(2).to_string())
rr = [stats.spearmanr(x.volts, x.mean_last100)[0] for _, x in g.groupby(["neuron", "carrier"]) if x.volts.nunique() >= 3]
print("Spearman rho of the last-100-ms mean dF on voltage: %d neuron-frequencies, median %.2f, positive in %d" % (len(rr), np.nanmedian(rr), np.sum(np.array(rr) > 0)))
neuron carrier volts n_ap n_ap_last100 mean_last100 pre_sd mean_post50
212 001032/10_1 2000.0 8.0 17 6 4.96 0.84 0.89
392 001032/11_1 20000.0 18.0 6 3 0.64 0.89 -1.95
452 001032/11_1 2000.0 18.0 3 1 0.33 0.97 1.38
500 001032/11_1 2000.0 18.0 3 1 0.59 0.88 -1.00
512 001032/11_1 2000.0 18.0 4 1 1.26 1.03 -1.13
572 001032/12_1 2000.0 15.0 0 0 9.82 0.89 7.87
656 001032/13_1 20000.0 14.0 10 2 -0.04 0.89 -1.87
692 001032/13_1 20000.0 14.0 4 1 1.05 0.91 -2.11
776 001032/13_1 2000.0 6.0 0 0 -0.48 0.89 -0.85
860 001032/14_1 20000.0 14.0 10 4 0.77 0.90 -2.43
956 001032/14_1 2000.0 4.0 4 2 -0.12 1.05 -1.03
980 001032/2_1 2000.0 10.0 8 3 3.27 1.00 NaN
992 001032/14_1 2000.0 4.0 7 3 1.63 0.98 -0.46
1040 001032/2_1 2000.0 10.0 5 2 5.01 0.89 NaN
1340 001032/3_1 2000.0 7.5 7 2 3.00 1.05 NaN
1532 001032/5_1 2000.0 8.0 8 3 2.53 0.94 0.01
1544 001032/5_1 2000.0 8.0 7 2 1.90 0.84 NaN
1712 001032/6_1 2000.0 5.0 2 1 0.87 1.04 -0.67
1928 001032/7_1 2000.0 8.0 13 4 6.02 0.94 -0.70
2096 001032/9_1 2000.0 7.5 8 3 4.81 0.97 4.06
2360 001033/10_2 20000.0 17.0 10 3 0.05 0.83 -3.54
2552 001033/10_2 2000.0 4.8 3 1 0.36 0.83 -1.21
2708 001033/2_2 2000.0 7.0 7 2 2.13 0.91 -0.12
2876 001033/3_2 2000.0 10.0 2 1 7.43 0.76 2.65
3548 001033/6_2 20000.0 20.0 17 8 3.63 0.76 -0.13
3716 001033/6_2 2000.0 9.0 7 3 2.35 0.79 -0.38
3848 001033/7_2 20000.0 16.0 0 0 -1.55 0.83 -1.73
3968 001033/7_2 2000.0 5.4 7 3 -0.05 0.71 -0.72
4040 001033/8_2 2000.0 10.0 1 0 2.67 0.96 -0.08
4220 001033/9_2 20000.0 12.0 7 3 1.18 0.99 -1.82
4388 001033/9_2 2000.0 4.0 10 3 2.20 1.07 -1.13
4544 001035/1_4 2000.0 6.4 5 2 2.28 0.71 -0.63
4712 001035/10_4 20000.0 18.0 0 0 -6.14 0.56 -11.04
4868 001035/11_4 2000.0 5.0 1 0 -0.05 0.77 -1.30
4952 001035/12_4 20000.0 10.0 0 0 -0.75 0.29 -1.23
4964 001035/12_4 20000.0 10.0 0 0 -0.72 1.34 -0.50
5120 001035/13_4 20000.0 10.0 4 3 0.69 0.96 -0.50
5216 001035/13_4 2000.0 4.0 2 0 0.42 0.88 -0.51
5300 001035/18_4 20000.0 8.8 3 0 -0.36 0.83 -0.47
5456 001035/2_4 20000.0 19.0 8 2 -0.60 0.75 -3.51
5528 001035/2_4 2000.0 4.5 10 4 1.49 0.89 -0.21
6128 001035/3_4 20000.0 20.0 8 1 2.00 0.77 -3.17
6224 001035/3_4 2000.0 6.0 8 3 4.71 0.75 -0.20
6296 001035/4_4 2000.0 4.6 3 1 2.00 0.77 0.47
6428 001035/5_4 2000.0 4.8 10 4 1.37 1.04 -1.05
6524 001035/6_4 20000.0 12.0 1 1 0.06 0.77 -0.32
6788 001035/6_4 2000.0 9.0 3 1 4.46 0.91 0.95
7052 001035/7_4 20000.0 20.0 8 0 0.18 0.78 -4.80
7172 001035/8_4 20000.0 10.0 10 5 0.60 0.82 -1.87
7256 001035/8_4 2000.0 3.8 5 2 0.11 0.80 -1.20
7364 001035/9_4 20000.0 11.0 3 2 0.39 0.87 -1.14
7484 001035/9_4 2000.0 3.0 6 2 1.34 0.84 -0.38
sessions meeting O4 (a) at any voltage: 2 in 2 of 32 neurons
neuron carrier volts n_ap mean_last100 pre_sd t_end
572 001032/12_1 2000.0 15.0 0 9.82 0.89 388.48
3044 001033/3_2 2000.0 9.0 5 6.74 0.75 388.48
Spearman rho of the last-100-ms mean dF on voltage: 45 neuron-frequencies, median 0.93, positive in 42
Result: not reproduced. Of the 12 to 14 neuron-frequencies tested above the voltage of their highest firing rate under the crop ROI, none shows sustained depolarization (O4, either level) at its highest voltage. Under ROI (a) 2 of 12 to 15 do, one of which falls below the O4 level within 50 ms of the stimulus end; under ROI (b) 1 of 4 to 13 under O4 (b) and none under O4 (a). Most neurons still fire in the last 100 ms at their highest voltage. Two sessions at any voltage meet O4 (a) under the crop ROI (001032/12_1 at 15 V, with no AP at all and a mean dF of 9.82 % in the last 100 ms, and 001033/3_2 at 9 V), and the mean dF in the last 100 ms rises with voltage (Spearman rho median 0.93, positive in 42 of 45 neuron-frequencies), so a depolarization builds up, but the block the claim describes is rare. 001032/12_1 stays depolarized after the stimulus (7.87 % over the next 50 ms), against the claim's immediate reversal. The traces are not corrected for bleaching (F2); the mean dF in the 50 ms after the stimulus is mostly below 0, so bleaching lowers the late dF by up to a few percent, which works against finding the block.
Thresholds (fig1-20khz-needs-twice-the-field, fig2d-5ap-above-1ap, fig2d-modulated-thresholds-higher)¶
Per neuron and stimulus (2 or 20 kHz, modulated or not), the lowest labelled voltage with at least 1 and at least 5 APs in the 250-ms stimulus (F5). Neurons that never reached the count are left out of that mean (counted below). Ratios do not depend on the conversion of volts to V/cm, but they depend on the files' labelled voltages being in the same proportion to the true amplitude at 2 and 20 kHz, which the files say they are not ("calibration factors are frequency-dependent").
def thresholds(g):
return g.groupby(["neuron", "stim"]).apply(lambda x: pd.Series(dict(t1=ts.threshold(x.volts, x.n_ap, 1), t5=ts.threshold(x.volts, x.n_ap, 5),
vmax=x.volts.max())), include_groups=False).reset_index()
TH = {r: thresholds(rd(FL, r)) for r in RD}
th = TH[("c", "a", "a")]
print("ROI c, F0 a, O3 a: neurons tested, reaching 1 AP and 5 APs, mean thresholds (labelled V)")
print(th.groupby("stim").agg(tested=("neuron", "size"), n1=("t1", "count"), t1_mean=("t1", "mean"), t1_sem=("t1", "sem"),
n5=("t5", "count"), t5_mean=("t5", "mean"), t5_sem=("t5", "sem")).round(2).to_string())
ROI c, F0 a, O3 a: neurons tested, reaching 1 AP and 5 APs, mean thresholds (labelled V)
tested n1 t1_mean t1_sem n5 t5_mean t5_sem
stim
mod_20k 3 3 14.00 2.08 2 17.20 2.80
mod_2k 14 13 5.38 0.45 8 6.15 0.63
pure_20k 17 15 9.85 0.70 11 13.55 0.85
pure_2k 28 27 4.22 0.34 22 5.00 0.36
fig1-20khz-needs-twice-the-field¶
All excitation thresholds at 20 kHz were about 2-fold higher than at 2 kHz.
Ratio of the mean 20-kHz to the mean 2-kHz threshold for each stimulus type and AP count (expected 2 +- 0.2), and two-tailed unpaired t tests: (a) one per type and count, each p < 0.0001, (b) one pooling the four.
rows = []
for r in RD:
th = TH[r]
d = dict(o1=r[0], o2=r[1], o3=r[2])
a20, a2, ps = [], [], []
for typ in ["pure", "mod"]:
for k in ["t1", "t5"]:
x20 = th[th.stim == f"{typ}_20k"][k].dropna(); x2 = th[th.stim == f"{typ}_2k"][k].dropna()
d[f"ratio_{typ}_{k}"] = x20.mean() / x2.mean()
p = stats.ttest_ind(x20, x2).pvalue if len(x20) > 1 and len(x2) > 1 else np.nan
d[f"p_{typ}_{k}"] = p; d[f"n_{typ}_{k}"] = f"{len(x20)}/{len(x2)}"
a20 += list(x20); a2 += list(x2); ps.append(p)
d["p_pooled"] = stats.ttest_ind(a20, a2).pvalue
ratios = [d[f"ratio_{t}_{k}"] for t in ["pure", "mod"] for k in ["t1", "t5"]]
d["ratios_ok"] = all(abs(x - 2) <= 0.2 for x in ratios if np.isfinite(x)) and all(np.isfinite(ratios))
d["holds_a"] = d["ratios_ok"] and all(np.isfinite(p) and p < 1e-4 for p in ps)
d["holds_b"] = d["ratios_ok"] and d["p_pooled"] < 1e-4
rows.append(d)
R20 = pd.DataFrame(rows)
with pd.option_context("display.float_format", "{:.3g}".format):
print(R20.to_string())
o1 o2 o3 ratio_pure_t1 p_pure_t1 n_pure_t1 ratio_pure_t5 p_pure_t5 n_pure_t5 ratio_mod_t1 p_mod_t1 n_mod_t1 ratio_mod_t5 p_mod_t5 n_mod_t5 p_pooled ratios_ok holds_a holds_b 0 c a a 2.33 5.12e-10 15/27 2.71 4.05e-12 11/22 2.6 1.12e-05 3/13 2.8 0.000198 2/8 7.85e-24 False False False 1 c a b 2.45 8.24e-10 16/27 2.71 4.05e-12 11/22 2.6 1.12e-05 3/13 2.8 0.000198 2/8 2.57e-24 False False False 2 c b a 2.33 5.12e-10 15/27 2.74 1.91e-12 11/22 2.6 1.12e-05 3/13 2.34 NaN 1/8 5.83e-24 False False False 3 c b b 2.45 8.24e-10 16/27 2.74 1.91e-12 11/22 2.6 1.12e-05 3/13 2.34 NaN 1/8 2.15e-24 False False False 4 a a a 2.48 2.96e-09 10/24 3.3 3.16e-07 7/10 2.37 4.71e-05 2/11 NaN NaN 0/1 4.6e-18 False False False 5 a a b 2.44 4.84e-07 10/20 3.52 2e-07 5/10 2.33 0.000704 2/8 NaN NaN 0/1 2.29e-14 False False False 6 a b a 2.46 4.97e-09 10/24 3 2.86e-06 6/11 2.22 0.000116 2/12 NaN NaN 0/1 2.5e-16 False False False 7 a b b 2.46 1.11e-06 10/19 3.23 5.03e-07 5/11 2.3 0.000699 2/8 NaN NaN 0/1 1.86e-13 False False False 8 b a a 2.55 4.8e-10 9/22 2.58 NaN 1/4 2.26 0.000368 2/9 NaN NaN 0/0 1.19e-15 False False False 9 b a b 2.35 0.000115 4/12 2.53 NaN 1/4 2.04 NaN 1/7 NaN NaN 0/0 1.91e-08 False False False 10 b b a 2.55 4.8e-10 9/22 2.14 NaN 1/5 2.15 0.00055 2/10 NaN NaN 0/0 1.26e-14 False False False 11 b b b 2.32 0.000189 4/11 2.11 NaN 1/5 NaN NaN 0/5 NaN NaN 0/0 2.51e-06 False False False
# Within neurons tested at both frequencies (unmodulated, 1 AP), ROI c, F0 a, O3 a: the ratio per neuron.
th = TH[("c", "a", "a")].pivot(index="neuron", columns="stim", values="t1")
both = th[["pure_2k", "pure_20k"]].dropna()
rr = both.pure_20k / both.pure_2k
print(f"{len(both)} neurons with both: ratio median {rr.median():.2f}, range {rr.min():.2f} to {rr.max():.2f}; paired t p = {stats.ttest_rel(both.pure_20k, both.pure_2k).pvalue:.2g}")
11 neurons with both: ratio median 3.36, range 1.79 to 4.17; paired t p = 1e-06
Result: blocked, metadata insufficient. The claim is a ratio of fields in V/cm at 20 and 2 kHz. The files give labelled generator voltages and say the delivered amplitude differs from the label by frequency-dependent calibration factors "reported separately", which are in none of the dandisets, so the ratio cannot be formed. Shown, not deciding: in labelled volts the 20-kHz thresholds are higher in every reading (pooled p 2.5e-06 to 2.2e-24; under the crop ROI unmodulated p 5.1e-10 to 8.2e-10 for 1 AP and 1.9e-12 to 4.1e-12 for 5 APs), with ratios of the means of 2.33 to 2.80 under the crop ROI (2.04 to 3.52 over all readings), outside 2 +- 0.2; within the 11 neurons tested at both frequencies the median ratio is 3.36. The 20-kHz modulated group has 3 neurons. The first round recorded this claim as partially reproduced.
fig2d-5ap-above-1ap¶
All thresholds to elicit 5 APs were higher than those for 1 AP.
rows = []
for r in RD:
th = TH[r]; d = dict(o1=r[0], o2=r[1], o3=r[2])
for s in ["pure_2k", "mod_2k", "pure_20k", "mod_20k"]:
x = th[th.stim == s]; d[s] = f"{x.t1.mean():.2f} / {x.t5.mean():.2f} (n {x.t1.count()}/{x.t5.count()})"
d[f"ok_{s}"] = x.t5.mean() > x.t1.mean()
d["holds"] = all(d[f"ok_{s}"] for s in ["pure_2k", "mod_2k", "pure_20k", "mod_20k"])
rows.append(d)
F5 = pd.DataFrame(rows); print(F5.to_string())
o1 o2 o3 pure_2k ok_pure_2k mod_2k ok_mod_2k pure_20k ok_pure_20k mod_20k ok_mod_20k holds 0 c a a 4.22 / 5.00 (n 27/22) True 5.38 / 6.15 (n 13/8) True 9.85 / 13.55 (n 15/11) True 14.00 / 17.20 (n 3/2) True True 1 c a b 4.22 / 5.00 (n 27/22) True 5.38 / 6.15 (n 13/8) True 10.36 / 13.55 (n 16/11) True 14.00 / 17.20 (n 3/2) True True 2 c b a 4.22 / 5.00 (n 27/22) True 5.38 / 6.15 (n 13/8) True 9.85 / 13.73 (n 15/11) True 14.00 / 14.40 (n 3/1) True True 3 c b b 4.22 / 5.00 (n 27/22) True 5.38 / 6.15 (n 13/8) True 10.36 / 13.73 (n 16/11) True 14.00 / 14.40 (n 3/1) True True 4 a a a 4.32 / 4.67 (n 24/10) True 5.07 / 8.00 (n 11/1) True 10.73 / 15.43 (n 10/7) True 12.00 / nan (n 2/0) False False 5 a a b 4.56 / 4.71 (n 20/10) True 5.15 / 8.00 (n 8/1) True 11.13 / 16.60 (n 10/5) True 12.00 / nan (n 2/0) False False 6 a b a 4.32 / 5.11 (n 24/11) True 5.40 / 8.00 (n 12/1) True 10.63 / 15.33 (n 10/6) True 12.00 / nan (n 2/0) False False 7 a b b 4.46 / 5.15 (n 19/11) True 5.21 / 8.00 (n 8/1) True 10.99 / 16.60 (n 10/5) True 12.00 / nan (n 2/0) False False 8 b a a 4.85 / 4.65 (n 22/4) False 5.31 / nan (n 9/0) False 12.33 / 12.00 (n 9/1) False 12.00 / nan (n 2/0) False False 9 b a b 4.74 / 4.75 (n 12/4) True 5.40 / nan (n 7/0) False 11.15 / 12.00 (n 4/1) True 11.00 / nan (n 1/0) False False 10 b b a 4.85 / 5.62 (n 22/5) True 5.59 / nan (n 10/0) False 12.33 / 12.00 (n 9/1) False 12.00 / nan (n 2/0) False False 11 b b b 4.81 / 5.70 (n 11/5) True 6.10 / nan (n 5/0) False 11.15 / 12.00 (n 4/1) True nan / nan (n 0/0) False False
Result: reproduced, resting on the crop ROI. Under the crop ROI the mean 5-AP threshold is above the mean 1-AP threshold for all four stimuli in all four F0 and O3 readings (F0 a, O3 a: 5.00 against 4.22 V for pure 2 kHz, 6.15 against 5.38 for modulated 2 kHz, 13.55 against 9.85 for pure 20 kHz and 17.20 against 14.00 for modulated 20 kHz, where 2 of 3 neurons reach 5 APs). Under the Otsu ROIs at most 1 neuron reaches 5 APs with modulated waves, so those comparisons cannot be made, and under ROI (b) with O3 (a) the 5-AP mean is below the 1-AP mean for pure 2 kHz (4.65 against 4.85 V, over 4 and 22 neurons) and pure 20 kHz (12.00 against 12.33 V, over 1 and 9), means over different neurons. Within a neuron the 5-AP threshold cannot be below its 1-AP threshold.
fig2d-modulated-thresholds-higher¶
Thresholds for sine waves modulated at 20 Hz were significantly higher than for pure sine waves at the same carrier frequency for 5 APs at 2 kHz and for 1 and 5 APs at 20 kHz, so modulation did not lower excitation thresholds; the 2-kHz unmodulated AP threshold was 57.7 +- 3.6 V/cm.
Two-tailed unpaired t tests (Figure 2D legend), with a paired test over neurons tested with both as a variant. The 2-kHz unmodulated threshold in V/cm from the labelled peak-to-peak volts under the two field readings (36 V/cm per volt of amplitude, half the peak-to-peak value; or per volt peak to peak), for the 1-AP and the 5-AP threshold.
rows = []
for r in RD:
th = TH[r]; d = dict(o1=r[0], o2=r[1], o3=r[2])
for car in ["2k", "20k"]:
for k in ["t1", "t5"]:
m = th[th.stim == f"mod_{car}"][k].dropna(); p_ = th[th.stim == f"pure_{car}"][k].dropna()
t = stats.ttest_ind(m, p_) if len(m) > 1 and len(p_) > 1 else None
d[f"{car}_{k}"] = f"{m.mean():.2f} vs {p_.mean():.2f}, p={t.pvalue:.3g}" if t else f"n={len(m)}/{len(p_)}"
d[f"up_{car}_{k}"] = bool(t is not None and m.mean() > p_.mean() and t.pvalue < 0.05)
d[f"down_{car}_{k}"] = bool(t is not None and m.mean() < p_.mean() and t.pvalue < 0.05)
x = th[th.stim == "pure_2k"]
for fr, fac in ts.FIELD_READINGS.items():
for k in ["t1", "t5"]:
d[f"Vcm_{fr}_{k}"] = x[k].mean() * fac * ts.V_TO_VCM
d["value_ok"] = any(abs(d[f"Vcm_{fr}_{k}"] - 57.7) <= 3.6 for fr in ts.FIELD_READINGS for k in ["t1", "t5"])
d["tests_ok"] = d["up_2k_t5"] and d["up_20k_t1"] and d["up_20k_t5"] and not d["down_2k_t1"]
rows.append(d)
MOD = pd.DataFrame(rows)
print(MOD[["o1", "o2", "o3", "2k_t1", "2k_t5", "20k_t1", "20k_t5", "tests_ok"]].to_string())
print(MOD[["o1", "o2", "o3"] + [c for c in MOD if c.startswith("Vcm")] + ["value_ok"]].round(1).to_string())
o1 o2 o3 2k_t1 2k_t5 20k_t1 20k_t5 tests_ok 0 c a a 5.38 vs 4.22, p=0.0564 6.15 vs 5.00, p=0.118 14.00 vs 9.85, p=0.0337 17.20 vs 13.55, p=0.133 False 1 c a b 5.38 vs 4.22, p=0.0564 6.15 vs 5.00, p=0.118 14.00 vs 10.36, p=0.102 17.20 vs 13.55, p=0.133 False 2 c b a 5.38 vs 4.22, p=0.0564 6.15 vs 5.00, p=0.118 14.00 vs 9.85, p=0.0337 n=1/11 False 3 c b b 5.38 vs 4.22, p=0.0564 6.15 vs 5.00, p=0.118 14.00 vs 10.36, p=0.102 n=1/11 False 4 a a a 5.07 vs 4.32, p=0.236 n=1/10 12.00 vs 10.73, p=0.538 n=0/7 False 5 a a b 5.15 vs 4.56, p=0.467 n=1/10 12.00 vs 11.13, p=0.751 n=0/5 False 6 a b a 5.40 vs 4.32, p=0.0898 n=1/11 12.00 vs 10.63, p=0.514 n=0/6 False 7 a b b 5.21 vs 4.46, p=0.354 n=1/11 12.00 vs 10.99, p=0.72 n=0/5 False 8 b a a 5.31 vs 4.85, p=0.53 n=0/4 12.00 vs 12.33, p=0.855 n=0/1 False 9 b a b 5.40 vs 4.74, p=0.442 n=0/4 n=1/4 n=0/1 False 10 b b a 5.59 vs 4.85, p=0.306 n=0/5 12.00 vs 12.33, p=0.855 n=0/1 False 11 b b b 6.10 vs 4.81, p=0.176 n=0/5 n=0/4 n=0/1 False o1 o2 o3 Vcm_amplitude_t1 Vcm_amplitude_t5 Vcm_peak_to_peak_t1 Vcm_peak_to_peak_t5 value_ok 0 c a a 76.0 90.1 152.0 180.2 False 1 c a b 76.0 90.1 152.0 180.2 False 2 c b a 76.0 90.1 152.0 180.2 False 3 c b b 76.0 90.1 152.0 180.2 False 4 a a a 77.8 84.1 155.5 168.1 False 5 a a b 82.2 84.8 164.3 169.6 False 6 a b a 77.8 92.0 155.5 183.9 False 7 a b b 80.3 92.6 160.7 185.2 False 8 b a a 87.2 83.7 174.4 167.4 False 9 b a b 85.4 85.5 170.7 171.0 False 10 b b a 87.2 101.2 174.4 202.3 False 11 b b b 86.6 102.6 173.1 205.2 False
# Paired variant: neurons tested with both, ROI c, F0 a, O3 a.
th = TH[("c", "a", "a")]
for car in ["2k", "20k"]:
for k in ["t1", "t5"]:
pv = th.pivot(index="neuron", columns="stim", values=k)[[f"mod_{car}", f"pure_{car}"]].dropna()
if len(pv) > 1:
t = stats.ttest_rel(pv[f"mod_{car}"], pv[f"pure_{car}"])
print(f"{car} {k}: {len(pv)} neurons, modulated higher in {(pv[f'mod_{car}'] > pv[f'pure_{car}']).sum()}, "
f"lower in {(pv[f'mod_{car}'] < pv[f'pure_{car}']).sum()}, mean difference {(pv[f'mod_{car}'] - pv[f'pure_{car}']).mean():.2f} V, paired p = {t.pvalue:.3g}")
else:
print(f"{car} {k}: {len(pv)} neurons with both")
# Which neurons were tested with 20-kHz modulated waves, and how far up.
print(th[th.stim == "mod_20k"].to_string())
2k t1: 12 neurons, modulated higher in 11, lower in 0, mean difference 0.64 V, paired p = 7.94e-05
2k t5: 8 neurons, modulated higher in 4, lower in 0, mean difference 0.44 V, paired p = 0.0893
20k t1: 3 neurons, modulated higher in 3, lower in 0, mean difference 3.17 V, paired p = 0.113
20k t5: 2 neurons, modulated higher in 1, lower in 0, mean difference 2.00 V, paired p = 0.5
neuron stim t1 t5 vmax
23 001033/10_2 mod_20k 13.0 14.4 17.0
31 001033/6_2 mod_20k 18.0 20.0 20.0
38 001033/9_2 mod_20k 11.0 NaN 12.0
Result: partially reproduced. Under the crop ROI modulated thresholds are higher than unmodulated ones in every comparison, but only one of the paper's three significant differences is significant here: 20 kHz with 1 AP, p = 0.034 under O3 (a) (0.10 under O3 b), with 3 modulated neurons. 20 kHz with 5 APs gives p = 0.13 (2 neurons) and 2 kHz with 5 APs p = 0.12 (paper 0.0079). 2 kHz with 1 AP does not differ the other way (p = 0.056, modulated higher). Paired within neurons, the 2-kHz 1-AP threshold is higher with modulation in 11 of 12 (p = 7.9e-05). Under the Otsu ROIs no test reaches p < 0.05. The mean unmodulated 2-kHz threshold is 4.22 V peak to peak, 76.0 V/cm (1 AP) or 90.1 V/cm (5 APs) if 36 V/cm applies to the amplitude and twice that if to the peak-to-peak value, none within 57.7 +- 3.6; the files say the delivered amplitude is "somewhat smaller" than the label and the calibration is not deposited, so this part cannot be judged. With the first round's detector (AP detector variants) the 20-kHz tests give p = 0.021 and 0.031 and the 2-kHz 5-AP test p = 0.09.
fig2-modulated-locks-to-beats¶
With 20-Hz modulation, firing above threshold locked to the beats, one AP near the peak of each of the five beats and at most two per beat (ten in all), and modulated sine waves triggered no sustained firing or sustained depolarization.
The files hold no stimulus waveform, so the envelope phase takes the claims file's two readings: (a) beat peaks 25 ms into each 50-ms beat, (b) at 0, 50, 100, 150 and 200 ms. Within 12.5 ms of a peak, pooled over sessions with at least 5 APs (and per session as a variant). The phase locking itself, whatever the phase, is shown by the vector strength of AP times at 20 Hz with a Rayleigh test.
def beats(g, o4, base=None):
base = FL if base is None else base
g = g[g.kind == "mod"].copy()
medh = rd(base, (g.o1.iloc[0], g.o2.iloc[0], g.o3.iloc[0])).groupby("neuron").ap_h.apply(
lambda s: np.median(np.concatenate(s.values)) if len(np.concatenate(s.values)) else np.nan)
d = dict(sessions=len(g), max_aps=int(g.n_ap.max()), over10=int((g.n_ap > 10).sum()))
sus = [ts.sustained(row, medh[row.neuron], o4)[0] for _, row in g.iterrows()]
d["sustained"] = int(np.sum(sus))
s5 = g[g.n_ap >= 5]
ph = np.concatenate([(x - 50.0) % 50.0 for x in s5.ap_t]) if len(s5) else np.array([])
d["sessions_5ap"] = len(s5); d["aps"] = len(ph)
d["near_a"] = np.mean(np.abs(ph - 25) <= 12.5) if len(ph) else np.nan
d["near_b"] = np.mean((ph <= 12.5) | (ph >= 37.5)) if len(ph) else np.nan
per_a = [np.mean(np.abs(((x - 50) % 50) - 25) <= 12.5) > 0.5 for x in s5.ap_t]
d["sessions_major_a"] = int(np.sum(per_a))
ang = 2 * np.pi * ph / 50.0
vs = np.abs(np.mean(np.exp(1j * ang))) if len(ph) else np.nan
d["vector_strength"] = vs; d["rayleigh_p"] = np.exp(-len(ph) * vs ** 2) if len(ph) else np.nan
d["mean_phase_ms"] = (np.angle(np.mean(np.exp(1j * ang))) % (2 * np.pi)) / (2 * np.pi) * 50 if len(ph) else np.nan
d["holds_a"] = d["over10"] == 0 and d["near_a"] > 0.5 and d["sustained"] == 0
d["holds_b"] = d["over10"] == 0 and d["near_b"] > 0.5 and d["sustained"] == 0
return d
BT = pd.DataFrame([dict(o1=r[0], o2=r[1], o3=r[2], o4=o4, **beats(rd(FL, r), o4)) for r in RD for o4 in "ab"])
with pd.option_context("display.float_format", "{:.3g}".format):
print(BT.to_string())
o1 o2 o3 o4 sessions max_aps over10 sustained sessions_5ap aps near_a near_b sessions_major_a vector_strength rayleigh_p mean_phase_ms holds_a holds_b 0 c a a a 84 8 0 1 16 84 0.964 0.0357 16 0.735 2.03e-20 26.6 False False 1 c a a b 84 8 0 0 16 84 0.964 0.0357 16 0.735 2.03e-20 26.6 True False 2 c a b a 84 8 0 1 16 84 0.964 0.0357 16 0.735 2.03e-20 26.6 False False 3 c a b b 84 8 0 0 16 84 0.964 0.0357 16 0.735 2.03e-20 26.6 True False 4 c b a a 84 8 0 1 15 79 0.962 0.038 15 0.726 8.03e-19 26.3 False False 5 c b a b 84 8 0 0 15 79 0.962 0.038 15 0.726 8.03e-19 26.3 True False 6 c b b a 84 8 0 1 15 79 0.962 0.038 15 0.726 8.03e-19 26.3 False False 7 c b b b 84 8 0 0 15 79 0.962 0.038 15 0.726 8.03e-19 26.3 True False 8 a a a a 84 5 0 0 1 5 1 0 1 0.919 0.0147 28.8 True False 9 a a a b 84 5 0 0 1 5 1 0 1 0.919 0.0147 28.8 True False 10 a a b a 84 5 0 0 1 5 1 0 1 0.919 0.0147 28.8 True False 11 a a b b 84 5 0 0 1 5 1 0 1 0.919 0.0147 28.8 True False 12 a b a a 84 5 0 0 1 5 1 0 1 0.919 0.0147 28.8 True False 13 a b a b 84 5 0 0 1 5 1 0 1 0.919 0.0147 28.8 True False 14 a b b a 84 5 0 0 1 5 1 0 1 0.919 0.0147 28.8 True False 15 a b b b 84 5 0 0 1 5 1 0 1 0.919 0.0147 28.8 True False 16 b a a a 84 3 0 0 0 0 NaN NaN 0 NaN NaN NaN False False 17 b a a b 84 3 0 0 0 0 NaN NaN 0 NaN NaN NaN False False 18 b a b a 84 3 0 0 0 0 NaN NaN 0 NaN NaN NaN False False 19 b a b b 84 3 0 0 0 0 NaN NaN 0 NaN NaN NaN False False 20 b b a a 84 3 0 0 0 0 NaN NaN 0 NaN NaN NaN False False 21 b b a b 84 3 0 0 0 0 NaN NaN 0 NaN NaN NaN False False 22 b b b a 84 3 0 0 0 0 NaN NaN 0 NaN NaN NaN False False 23 b b b b 84 3 0 0 0 0 NaN NaN 0 NaN NaN NaN False False
# AP phase histogram within the 50-ms beat (ROI c, F0 a, O3 a), sessions with at least 5 APs, by carrier.
g = rd(FL, ("c", "a", "a")); g = g[(g.kind == "mod") & (g.n_ap >= 5)]
for f, x in g.groupby("carrier"):
ph = np.concatenate([(v - 50.0) % 50.0 for v in x.ap_t])
print(f"{int(f/1000)} kHz: {len(x)} sessions, {len(ph)} APs; counts per 5-ms bin of the beat:", np.histogram(ph, bins=np.arange(0, 51, 5))[0].tolist())
print("AP counts per modulated session:", rd(FL, ("c", "a", "a")).query("kind == 'mod'").n_ap.value_counts().sort_index().to_dict())
# The late window of the modulated sessions (ROI c, F0 a, O3 a): mean dF in the last 100 ms of the stimulus.
m = rd(FL, ("c", "a", "a")); medh_m = m.groupby("neuron").ap_h.apply(lambda s: np.median(np.concatenate(s.values)) if len(np.concatenate(s.values)) else np.nan)
m = m[m.kind == "mod"]
lv_b = m.neuron.map(0.5 * medh_m)
print("modulated sessions: %d; last-100-ms mean dF above 5 %%: %d; above the O4 (a) level: %d; above the O4 (b) level: %d; of these with no AP there: %d (a), %d (b)"
% (len(m), int((m.mean_last100 > 5).sum()), int((m.mean_last100 > 5 * m.pre_sd).sum()), int((m.mean_last100 > lv_b).sum()),
int(((m.mean_last100 > 5 * m.pre_sd) & (m.n_ap_last100 == 0)).sum()), int(((m.mean_last100 > lv_b) & (m.n_ap_last100 == 0)).sum())))
m = m.assign(o4a_level=5 * m.pre_sd, o4b_level=lv_b, median_ap_h=m.neuron.map(medh_m))
print(m.nlargest(5, "mean_last100")[["neuron", "carrier", "volts", "n_ap", "n_ap_last100", "mean_last100", "o4a_level", "o4b_level", "median_ap_h", "t_end"]].round(2).to_string())
2 kHz: 13 sessions, 69 APs; counts per 5-ms bin of the beat: [0, 0, 2, 14, 15, 18, 15, 3, 2, 0] 20 kHz: 3 sessions, 15 APs; counts per 5-ms bin of the beat: [0, 0, 0, 0, 3, 7, 5, 0, 0, 0]
AP counts per modulated session: {0: 23, 1: 11, 2: 10, 3: 17, 4: 7, 5: 14, 6: 1, 8: 1}
modulated sessions: 84; last-100-ms mean dF above 5 %: 1; above the O4 (a) level: 2; above the O4 (b) level: 0; of these with no AP there: 1 (a), 0 (b)
neuron carrier volts n_ap n_ap_last100 mean_last100 o4a_level o4b_level median_ap_h t_end
1124 001032/2_1 2000.0 11.0 3 0 5.96 4.72 6.10 12.20 291.28
308 001032/10_1 2000.0 8.0 8 4 4.52 4.43 6.17 12.34 388.48
296 001032/10_1 2000.0 5.0 2 1 4.29 4.46 6.17 12.34 388.48
248 001032/10_1 2000.0 5.0 2 1 2.58 4.12 6.17 12.34 388.48
1964 001032/7_1 2000.0 5.5 4 1 1.98 4.19 7.22 14.45 388.48
Result: reproduced, resting on the envelope phase and the O4 level. No modulated session has more than 10 APs (at most 8, of 84 sessions). In the 16 sessions with at least 5 APs (crop ROI), 96.4 % of APs fall within 12.5 ms of beat peaks 25 ms into each beat (envelope reading a) and 3.6 % near peaks at the start of each beat (reading b); the APs lock to 26.6 ms of the 50-ms beat (vector strength 0.735, Rayleigh p = 2.0e-20), and the phase histograms peak 15 to 35 ms into the beat at both carriers. Reading (a), the wave starting at minimal amplitude, is adopted: the paper's Methods say the modulated strobe waves start at minimal amplitude, and under reading (b) the APs would sit at the envelope's minima, where the field is near zero. No modulated session shows sustained depolarization under O4 (b); under O4 (a) one does under the crop ROI (001032/2_1 at 11 V, 5.96 % with no AP in its last 100 ms, a 900-frame session cut at 291 ms), none under the Otsu ROIs. It is also the only modulated session with a last-100-ms mean dF above 5 %. Under the Otsu ROIs at most 1 session reaches 5 APs, so the locking is shown under the crop ROI.
fig2c-pure-more-aps¶
In the same neuron, the threshold was typically lower for unmodulated sine waves, and at the same field strength pure sine waves triggered more APs than modulated ones.
Pairs are a neuron's modulated and unmodulated sessions at the same carrier and labelled voltage (AP counts averaged over repeated sessions); 1-AP thresholds compared within neurons tested with both.
rows = []
for r in RD:
g = rd(FL, r)
c = g.groupby(["neuron", "carrier", "volts", "kind"]).n_ap.mean().unstack("kind").dropna()
th = TH[r]; out = dict(o1=r[0], o2=r[1], o3=r[2], pairs=len(c), pure_more=int((c.sine > c["mod"]).sum()),
mod_more=int((c["mod"] > c.sine).sum()), equal=int((c["mod"] == c.sine).sum()),
paired_t_p=stats.ttest_rel(c.sine, c["mod"]).pvalue if len(c) > 1 else np.nan,
mean_diff=(c.sine - c["mod"]).mean())
lo = hi = n = 0
for car in ["2k", "20k"]:
pv = th.pivot(index="neuron", columns="stim", values="t1")
if f"mod_{car}" in pv and f"pure_{car}" in pv:
pv = pv[[f"mod_{car}", f"pure_{car}"]].dropna(); n += len(pv)
lo += int((pv[f"pure_{car}"] < pv[f"mod_{car}"]).sum()); hi += int((pv[f"pure_{car}"] > pv[f"mod_{car}"]).sum())
out.update(neurons_both=n, pure_threshold_lower=lo, pure_threshold_higher=hi)
out["holds"] = out["pure_more"] > out["mod_more"] and lo > hi
rows.append(out)
PC = pd.DataFrame(rows)
with pd.option_context("display.float_format", "{:.3g}".format):
print(PC.to_string())
o1 o2 o3 pairs pure_more mod_more equal paired_t_p mean_diff neurons_both pure_threshold_lower pure_threshold_higher holds 0 c a a 52 44 2 6 1.24e-08 3.35 15 14 0 True 1 c a b 52 44 2 6 1.31e-08 3.36 15 14 0 True 2 c b a 52 44 2 6 2.58e-08 3.31 15 14 0 True 3 c b b 52 44 2 6 2.41e-08 3.35 15 14 0 True 4 a a a 52 32 6 14 1.78e-05 1.34 11 9 0 True 5 a a b 52 24 2 26 7.24e-05 1.23 9 7 0 True 6 a b a 52 29 6 17 6.62e-05 1.26 11 9 0 True 7 a b b 52 25 1 26 4.98e-05 1.24 9 7 0 True 8 b a a 52 19 9 24 0.0353 0.441 10 5 3 True 9 b a b 52 10 8 34 0.635 0.0913 6 3 0 True 10 b b a 52 16 10 26 0.105 0.306 10 7 2 True 11 b b b 52 9 6 37 0.521 0.111 3 0 1 False
Result: reproduced, resting on the crop ROI. Under the crop ROI unmodulated waves gave more APs than modulated ones at the same carrier and voltage in 44 of 52 pairs against 2 (mean difference 3.35 APs, paired p = 1.2e-08), and the unmodulated 1-AP threshold was lower in 14 of 15 neurons tested with both (higher in none). Both parts hold in 11 of the 12 readings. Under ROI (b) with F0 (b) and O3 (b) the counts still favour pure waves (9 pairs against 6), but of the 3 neurons with both thresholds the pure one is lower in none and higher in 1.
fig3c-sustained-firing-2s¶
In the low-conductivity solution, a 2-s unmodulated 2-kHz sine wave ramped up over 250 ms evoked firing for its whole duration; the firing rate rose with field until replaced by a firing block with sustained depolarization, which developed before the end of the ramp.
001030: four neurons (c5 to c8) in NaGlu, 2-s waves starting at 25 ms with a 250-ms ramp (session descriptions). Bleaching is corrected "by subtracting a trace without APs" (Methods, Recording of the optical APs; F2). Only one session in 001030 has no AP (c7 at 1.0 V), so the corrections are: (i) c7's AP-free trace subtracted from every session (the only one there is), (ii) none, as a variant. Per neuron, the criterion's three parts: APs in the last 250 ms of the stimulus at some voltage; the AP count rising with voltage up to some voltage; and at the highest voltage no AP and mean dF above the O4 level from 275 ms (250 ms after onset) to the end. The figure shows one neuron, so the criterion is taken to hold when one neuron meets all three (a reading made here).
R30 = [r for r in REC["001030"] if "error" not in r]
ref = [r for r in R30 if ts.neuron_of_001030(r["path"]) == "c7" and abs(r["stim"]["volts"] - 1.0) < 1e-6][0]
def fig3c(o1, o2, o3, corr, width_ref="both"):
rows = []
for r in R30:
dF, info = r["traces"][(o1, o2)]
t = r["t"]
if corr == "c7_free":
if r is ref: continue
rd_, _ = ref["traces"][(o1, o2)]; n = min(len(dF), len(rd_)); dF = dF[:n] - rd_[:n]; t = t[:n]
pre_sd = float(dF[t < 25].std())
pt, ph, pw = ts.detect_aps(t, dF, 25, 2025, pre_sd, o3, width_ref)
late = (t >= 275) & (t < 2025)
rows.append(dict(neuron=ts.neuron_of_001030(r["path"]), volts=r["stim"]["volts"], n_ap=len(pt), last250=int((pt >= 1775).sum()),
aps_after_275=int((pt >= 275).sum()), mean_after_275=float(dF[late].mean()), pre_sd=pre_sd,
median_h=np.median(ph) if len(ph) else np.nan, t_end=float(t[-1])))
return pd.DataFrame(rows).sort_values(["neuron", "volts"])
def judge(df, o4):
out = {}
for n, x in df.groupby("neuron"):
x = x.sort_values("volts"); top = x.iloc[-1]
medh = np.nanmedian(x.median_h)
lvl = 5 * top.pre_sd if o4 == "a" else 0.5 * medh
rises = len(x) > 1 and x.n_ap.iloc[: int(np.argmax(x.n_ap.values)) + 1].is_monotonic_increasing and x.n_ap.values.argmax() > 0
out[n] = dict(late_aps=bool((x.last250 > 0).any()), rises=bool(rises), block=bool(top.aps_after_275 == 0 and top.mean_after_275 > lvl),
top_volts=top.volts, top_mean=round(top.mean_after_275, 2), level=round(lvl, 2))
out[n]["all3"] = out[n]["late_aps"] and out[n]["rises"] and out[n]["block"]
return out
rows = []
for o1 in "cab":
for o2 in "ab":
for o3 in "ab":
for corr in ["c7_free", "none"]:
df = fig3c(o1, o2, o3, corr)
for o4 in "ab":
j = judge(df, o4)
rows.append(dict(o1=o1, o2=o2, o3=o3, corr=corr, o4=o4, **{f"{n}": ("".join("Y" if v[k] else "-" for k in ["late_aps", "rises", "block"]))
for n, v in j.items()}, holds=any(v["all3"] for v in j.values())))
F3C = pd.DataFrame(rows)
print("per neuron: late APs / count rises / block at the highest voltage (Y = met)")
print(F3C.to_string())
per neuron: late APs / count rises / block at the highest voltage (Y = met) o1 o2 o3 corr o4 c5 c6 c7 c8 holds 0 c a a c7_free a YY- YY- YY- YY- False 1 c a a c7_free b YY- YY- YY- YY- False 2 c a a none a YY- YY- YY- YY- False 3 c a a none b YY- YY- YY- YY- False 4 c a b c7_free a YY- YY- YY- YY- False 5 c a b c7_free b YY- YY- YY- YY- False 6 c a b none a -Y- -Y- YY- YY- False 7 c a b none b -Y- -Y- YY- YY- False 8 c b a c7_free a YY- YY- YY- YY- False 9 c b a c7_free b YY- YY- YY- YY- False 10 c b a none a -Y- YY- YY- YY- False 11 c b a none b -Y- YY- YY- YY- False 12 c b b c7_free a YY- YY- YY- YY- False 13 c b b c7_free b YY- YY- YY- YY- False 14 c b b none a -Y- -Y- YY- YY- False 15 c b b none b -Y- -Y- YY- YY- False 16 a a a c7_free a YY- YY- YY- -Y- False 17 a a a c7_free b YY- YY- YY- -Y- False 18 a a a none a YY- YY- YY- -Y- False 19 a a a none b YY- YY- YY- -Y- False 20 a a b c7_free a YY- YY- YY- --- False 21 a a b c7_free b YY- YY- YY- --- False 22 a a b none a -Y- -Y- YY- --- False 23 a a b none b -Y- -Y- YY- --- False 24 a b a c7_free a YY- YY- YY- -Y- False 25 a b a c7_free b YY- YY- YY- -Y- False 26 a b a none a -Y- -Y- YY- -Y- False 27 a b a none b -Y- -Y- YY- -Y- False 28 a b b c7_free a YY- YY- YY- --- False 29 a b b c7_free b YY- YY- YY- --- False 30 a b b none a -Y- -Y- YY- --- False 31 a b b none b -Y- -Y- YY- --- False 32 b a a c7_free a -Y- YY- YY- -Y- False 33 b a a c7_free b -Y- YY- YY- -Y- False 34 b a a none a -Y- -Y- YY- -Y- False 35 b a a none b -Y- -Y- YY- -Y- False 36 b a b c7_free a --- --- YY- --- False 37 b a b c7_free b --- --- YY- --- False 38 b a b none a -Y- -Y- YY- --- False 39 b a b none b -Y- -Y- YY- --- False 40 b b a c7_free a -Y- YY- YY- -Y- False 41 b b a c7_free b -Y- YY- YY- -Y- False 42 b b a none a -Y- -Y- YY- -Y- False 43 b b a none b -Y- -Y- YY- -Y- False 44 b b b c7_free a --- --- YY- --- False 45 b b b c7_free b --- --- YY- --- False 46 b b b none a -Y- -Y- YY- --- False 47 b b b none b -Y- -Y- YY- --- False
df = fig3c("c", "a", "a", "c7_free")
print(df.round(2).to_string())
print("sessions shorter than the stimulus (end before 2025 ms):", df[df.t_end < 2025][["neuron", "volts", "t_end"]].values.tolist())
# The APs after 275 ms at the two highest-voltage sessions that come closest to a block, with their
# prominence (rise above the higher of the two neighbouring troughs) and width, and the count under
# the first round's detector (width at half the prominence) and under the rise-within-5-ms test the
# round-1 verifier asked for (dF rises by at least 5 % within 5 ms before the peak; not a claims-file reading).
from scipy import signal as _sg
for n_, v_ in [("c7", 2.3), ("c8", 2.0)]:
r = [r for r in R30 if ts.neuron_of_001030(r["path"]) == n_ and abs(r["stim"]["volts"] - v_) < 1e-6][0]
dF0, _ = r["traces"][("c", "a")]; rd_, _ = ref["traces"][("c", "a")]; m = min(len(dF0), len(rd_)); dF = (dF0[:m] - rd_[:m]).astype(float); t = r["t"][:m]
pre_sd = float(dF[t < 25].std())
for wr in ["both", "prominence"]:
pt, ph, pw = ts.detect_aps(t, dF, 25, 2025, pre_sd, "a", wr)
sel = pt >= 275
idx = np.searchsorted(t, pt[sel]); prom = _sg.peak_prominences(dF, idx)[0] if len(idx) else []
print(f"{n_} {v_} V, width {wr}: {int(sel.sum())} APs after 275 ms at", np.round(pt[sel], 0).tolist(), "height", np.round(ph[sel], 1).tolist(),
"prominence", np.round(prom, 1).tolist(), "width", np.round(pw[sel], 1).tolist())
k = int(round(5.0 / np.median(np.diff(t))))
pk, _ = _sg.find_peaks(dF, height=5.0, distance=k)
rise = [p for p in pk if t[p] >= 275 and t[p] < 2025 and dF[p] - dF[max(0, p - k):p].min() >= 5.0]
print(f" rise of 5 % within 5 ms: {len(rise)} peaks after 275 ms; mean dF after 275 ms {dF[(t >= 275) & (t < 2025)].mean():.2f} %")
fig, ax = plt.subplots(4, 1, figsize=(10, 9), sharex=True)
for i, n in enumerate(["c5", "c6", "c7", "c8"]):
for r in sorted([r for r in R30 if ts.neuron_of_001030(r["path"]) == n], key=lambda r: r["stim"]["volts"]):
dF = r["traces"][("c", "a")][0]; rd_ = ref["traces"][("c", "a")][0]; m = min(len(dF), len(rd_))
ax[i].plot(r["t"][:m], dF[:m] - rd_[:m] + 15 * (r["stim"]["volts"] - 1), lw=.4, label=f'{r["stim"]["volts"]} V')
ax[i].set_ylabel(f"{n} dF (%), offset"); ax[i].legend(fontsize=6, loc="upper right")
ax[-1].set_xlabel("ms"); plt.tight_layout(); plt.show()
neuron volts n_ap last250 aps_after_275 mean_after_275 pre_sd median_h t_end 2 c5 1.00 14 0 11 1.73 0.97 6.84 2232.68 3 c5 1.50 64 4 57 3.83 0.95 9.84 1913.68 4 c5 2.00 41 2 31 3.48 0.85 6.51 1913.68 1 c6 1.00 11 2 10 0.04 1.00 10.28 2232.68 0 c6 1.25 30 3 27 1.70 0.87 10.41 2232.68 8 c6 1.45 39 5 36 2.82 0.91 9.51 2232.68 5 c7 1.25 17 2 16 1.12 0.86 15.21 2232.68 7 c7 1.45 70 8 68 3.35 0.79 14.74 2232.68 6 c7 1.70 89 10 83 5.22 0.85 13.45 2232.68 13 c7 2.00 88 10 80 6.13 0.82 11.03 2232.68 11 c7 2.30 8 1 1 6.72 0.89 13.40 2232.68 10 c8 1.00 29 3 26 0.39 0.78 11.35 2232.68 12 c8 1.25 45 5 41 2.86 0.80 11.66 2232.68 9 c8 1.50 50 6 45 4.00 0.90 10.77 2232.68 15 c8 1.75 46 6 41 3.94 0.84 8.82 2232.68 14 c8 2.00 8 0 6 3.63 0.82 6.45 2232.68 sessions shorter than the stimulus (end before 2025 ms): [['c5', 1.5, 1913.6787109375], ['c5', 2.0, 1913.6787109375]] c7 2.3 V, width both: 1 APs after 275 ms at [1964.0] height [8.0] prominence [4.4] width [8.2] c7 2.3 V, width prominence: 11 APs after 275 ms at [302.0, 346.0, 360.0, 452.0, 678.0, 709.0, 1183.0, 1647.0, 1653.0, 1671.0, 1964.0] height [9.2, 9.7, 10.6, 9.5, 8.7, 9.6, 7.9, 8.8, 7.8, 8.2, 8.0] prominence [3.2, 4.4, 6.4, 3.9, 3.4, 4.6, 2.8, 5.0, 2.7, 4.2, 4.4] width [2.7, 2.7, 2.3, 2.9, 2.3, 2.5, 2.2, 2.2, 2.4, 4.0, 2.9] rise of 5 % within 5 ms: 1 peaks after 275 ms; mean dF after 275 ms 6.72 % c8 2.0 V, width both: 6 APs after 275 ms at [852.0, 924.0, 988.0, 1017.0, 1134.0, 1266.0] height [6.7, 5.4, 7.5, 5.4, 5.8, 6.1] prominence [5.9, 4.1, 7.2, 4.3, 4.5, 4.5] width [2.3, 3.0, 4.9, 3.9, 2.4, 3.0] c8 2.0 V, width prominence: 6 APs after 275 ms at [852.0, 924.0, 988.0, 1017.0, 1134.0, 1266.0] height [6.7, 5.4, 7.5, 5.4, 5.8, 6.1] prominence [5.9, 4.1, 7.2, 4.3, 4.5, 4.5] width [2.2, 2.0, 3.4, 3.0, 2.2, 2.3] rise of 5 % within 5 ms: 36 peaks after 275 ms; mean dF after 275 ms 3.63 %
Result: partially reproduced. In all four neurons of 001030, under the crop ROI with the bleaching correction, some voltage gives APs in the last 250 ms of the 2-s wave (firing for the whole duration) and the AP count rises with voltage up to a point. No neuron, under any reading or correction, shows the third part, a block at the highest voltage, though c7 at 2.3 V comes close: 7 APs during the ramp and one after 275 ms (at 1,964 ms, 8.0 % high, prominence 4.4 %), on a mean dF of 6.72 %. The first round's detector counted 11 APs there, with prominences of 2.7 to 6.4 %. c8 at 2.0 V keeps 6 APs after 275 ms (prominence 4.1 to 7.2 %) on a mean dF of 3.63 %. At the top voltage the count falls (c7 from 88 to 8, c8 from 46 to 8) while the mean dF stays raised. Two c5 sessions end at 1,914 ms, before the stimulus ends.
figs1-di8-same-as-fluovolt¶
With Di-8-ANEPPS in place of FluoVolt, in eight neurons, 2-kHz sine waves near threshold elicited one or a few delayed APs, and stronger fields shortened the delay and induced sustained firing ending in sustained depolarization, as with FluoVolt.
001272, traces inverted (Di-8-ANEPPS emission falls with depolarization; F1). O3 reading (a) is used only if the inverted APs reach 5% (claims file); their heights are shown first. Per neuron: one to three APs at the lowest voltage with APs; Spearman rho of first-AP latency on voltage below 0; sustained firing (an AP in the last 100 ms of the stimulus, a standard choice) or sustained depolarization (O4) at a voltage above the lowest with APs.
D8 = TAB["001272"]
x = rd(D8, ("c", "a", "b")); hh = np.concatenate(x.ap_h.values)
print("inverted Di-8 AP heights under O3 (b), ROI c: n %d, median %.1f %%, 90th percentile %.1f %%, share above 5%%: %.0f %%" % (len(hh), np.median(hh), np.percentile(hh, 90), 100 * np.mean(hh > 5)))
print("pre-onset SD median %.2f %% (5 SD = %.1f %%)" % (x.pre_sd.median(), 5 * x.pre_sd.median()))
rows = []
for r in RD:
g = rd(D8, r)
medh = g.groupby("neuron").ap_h.apply(lambda s: np.median(np.concatenate(s.values)) if len(np.concatenate(s.values)) else np.nan)
for n, y in g.groupby("neuron"):
a = y[y.n_ap > 0]
if len(a) == 0:
rows.append(dict(o1=r[0], o2=r[1], o3=r[2], neuron=n, sessions=len(y))); continue
v0 = a.volts.min(); low = a[a.volts == v0]
rho = stats.spearmanr(a.volts, a.latency)[0] if a.volts.nunique() > 1 else np.nan
hi = y[y.volts > v0]
sf = any((row.n_ap_last100 > 0) or ts.sustained(row, medh[n], "a")[0] or ts.sustained(row, medh[n], "b")[0] for _, row in hi.iterrows())
rows.append(dict(o1=r[0], o2=r[1], o3=r[2], neuron=n, sessions=len(y), v_lowest=v0, aps_at_lowest=low.n_ap.mean(),
one_to_three=bool(1 <= low.n_ap.median() <= 3), latency_at_lowest=low.latency.mean(), rho=rho,
sustained_above=sf))
DI = pd.DataFrame(rows)
with pd.option_context("display.float_format", "{:.3g}".format):
print(DI.to_string())
summ = DI.groupby(["o1", "o2", "o3"]).apply(lambda z: pd.Series(dict(neurons=len(z), one_to_three=int(z.one_to_three.fillna(False).sum()),
rho_neg=int((z.rho < 0).sum()), sustained=int(z.sustained_above.fillna(False).sum()))), include_groups=False)
summ["holds_on_these"] = (summ.one_to_three > summ.neurons / 2) & (summ.rho_neg > summ.neurons / 2) & (summ.sustained > summ.neurons / 2)
print(summ.to_string())
inverted Di-8 AP heights under O3 (b), ROI c: n 206, median 4.2 %, 90th percentile 5.1 %, share above 5%: 13 %
pre-onset SD median 0.55 % (5 SD = 2.8 %)
o1 o2 o3 neuron sessions v_lowest aps_at_lowest one_to_three latency_at_lowest rho sustained_above
0 c a a 001272/2_20240806 20 5 1 True 188 1 True
1 c a a 001272/3_20240806 29 2 1 True 64.9 -0.0445 True
2 c a a 001272/4_20240806 16 5.7 1 True 212 0.5 True
3 c a b 001272/2_20240806 20 4 2 True 5.73 0.0183 True
4 c a b 001272/3_20240806 29 2 1 True 63.4 -0.0311 True
5 c a b 001272/4_20240806 16 5.7 1.5 True 167 -0.687 True
6 c b a 001272/2_20240806 20 NaN NaN NaN NaN NaN NaN
7 c b a 001272/3_20240806 29 2.5 2.5 True 76 0.289 True
8 c b a 001272/4_20240806 16 5.7 1 True 212 0.5 True
9 c b b 001272/2_20240806 20 5 4.8 False 61.6 0.246 True
10 c b b 001272/3_20240806 29 2 1 True 63.4 -0.121 True
11 c b b 001272/4_20240806 16 5.7 1.5 True 167 -0.763 True
12 a a a 001272/2_20240806 20 5 1.67 True 101 0.0567 True
13 a a a 001272/3_20240806 29 2 1 True 63.3 -0.212 True
14 a a a 001272/4_20240806 16 5.7 1 True 211 -0.0736 True
15 a a b 001272/2_20240806 20 NaN NaN NaN NaN NaN NaN
16 a a b 001272/3_20240806 29 2 1 True 63.3 -0.319 True
17 a a b 001272/4_20240806 16 5.7 1 True 212 0.264 True
18 a b a 001272/2_20240806 20 5 1 True 172 -0.738 True
19 a b a 001272/3_20240806 29 2 1 True 66.6 -0.0655 True
20 a b a 001272/4_20240806 16 5.7 1 True 210 0.0522 True
21 a b b 001272/2_20240806 20 NaN NaN NaN NaN NaN NaN
22 a b b 001272/3_20240806 29 2 1 True 66.6 -0.28 True
23 a b b 001272/4_20240806 16 6.5 1 True 233 NaN False
24 b a a 001272/2_20240806 20 4 1 True 208 -0.5 True
25 b a a 001272/3_20240806 29 2 1 True 72.5 0.277 True
26 b a a 001272/4_20240806 16 5.7 1 True 211 -0.0736 True
27 b a b 001272/2_20240806 20 NaN NaN NaN NaN NaN NaN
28 b a b 001272/3_20240806 29 NaN NaN NaN NaN NaN NaN
29 b a b 001272/4_20240806 16 5.7 1 True 212 0.264 True
30 b b a 001272/2_20240806 20 4 1 True 208 -1 True
31 b b a 001272/3_20240806 29 2 1 True 72.5 0.104 True
32 b b a 001272/4_20240806 16 5.7 1 True 210 0.0522 True
33 b b b 001272/2_20240806 20 NaN NaN NaN NaN NaN NaN
34 b b b 001272/3_20240806 29 NaN NaN NaN NaN NaN NaN
35 b b b 001272/4_20240806 16 6.5 1 True 233 NaN False
neurons one_to_three rho_neg sustained holds_on_these
o1 o2 o3
a a a 3 3 2 3 True
b 3 2 1 2 False
b a 3 3 2 3 True
b 3 2 1 1 False
b a a 3 3 2 3 True
b 3 1 0 1 False
b a 3 3 1 3 False
b 3 1 0 0 False
c a a 3 3 1 3 False
b 3 3 2 3 True
b a 3 2 0 2 False
b 3 2 2 3 True
fig, ax = plt.subplots(3, 1, figsize=(10, 7), sharex=True)
for i, n in enumerate(sorted(set(r["subject_id"] for r in REC["001272"] if "error" not in r))):
rr = sorted([r for r in REC["001272"] if "error" not in r and r["subject_id"] == n], key=lambda r: r["start"])
for k, r in enumerate(rr[::3]):
ax[i].plot(r["t"], r["traces"][("c", "a")][0] + 6 * k, lw=.4)
ax[i].axvline(50, c="k", lw=.5); ax[i].axvline(300, c="k", lw=.5); ax[i].set_ylabel(f"{n}\n-dF (%), offset")
ax[-1].set_xlabel("ms"); plt.tight_layout(); plt.show()
Result: partially reproduced. 001272 holds 3 Di-8-ANEPPS neurons, not the paper's 8. The inverted APs are small (median 4.2 % dF, 13 % of them above 5 %), so O3 reading (b) alone applies, as the claims file says. Under the crop ROI with O3 (b), most of the 3 neurons show one to three APs at their lowest voltage with APs (3 with F0 a, 2 with F0 b), a latency that falls with voltage (2 of 3, rho -0.03 and -0.69; -0.12 and -0.76 with F0 b) and sustained firing or depolarization at higher voltages (3 of 3). Under the Otsu ROIs with O3 (b) the latency falls in at most 1 of 3 and one or two neurons have no AP. The behaviour holds in the deposited neurons under the crop ROI (an assumption); the count of eight does not, five neurons not being deposited.
figs2-ttx-blocks-depolarization¶
1 uM tetrodotoxin abolished both the APs and the sustained depolarization evoked by intense unmodulated 2-kHz sine waves, while the vehicle did not.
- In TTX neurons the sessions before are those labelled vehicle and the sessions after those labelled TTX. Vehicle neurons have no label for the change of solution; their sessions fall into blocks separated by pauses of about 300 s, and the split is taken at the longest pause (an assumption; the pause is printed). Responses are compared at the voltages tested both before and after. A response is an AP in the stimulus or sustained depolarization (O4).
TT = TAB["001271"]
def ttx(r, o4, table=None):
g = rd(TT if table is None else table, r); g = g[g.kind == "sine"]
medh = g.groupby("neuron").ap_h.apply(lambda s: np.median(np.concatenate(s.values)) if len(np.concatenate(s.values)) else np.nan)
rows = []
for n, x in g.groupby("neuron"):
b, a, gap = ts.split_before_after(x)
common = sorted(set(b.volts) & set(a.volts))
if not common:
rows.append(dict(neuron=n, group="TTX" if x.ttx.any() else "vehicle", gap_s=gap, common_volts=0)); continue
def resp(y):
return [bool(row.n_ap > 0 or ts.sustained(row, medh[n], o4)[0]) for _, row in y[y.volts.isin(common)].iterrows()]
def susd(y):
return int(sum(ts.sustained(row, medh[n], o4)[0] for _, row in y[y.volts.isin(common)].iterrows()))
rb, ra = resp(b), resp(a)
rows.append(dict(neuron=n, group="TTX" if x.ttx.any() else "vehicle", gap_s=gap, common_volts=len(common),
before_resp=f"{sum(rb)}/{len(rb)}", after_resp=f"{sum(ra)}/{len(ra)}",
aps_before=int(b[b.volts.isin(common)].n_ap.sum()), aps_after=int(a[a.volts.isin(common)].n_ap.sum()),
sus_before=susd(b), sus_after=susd(a),
ok=(all(rb) and not any(ra)) if x.ttx.any() else all(ra),
ok_any=(any(rb) and not any(ra)) if x.ttx.any() else any(ra)))
return pd.DataFrame(rows)
rows = []
for r in RD:
for o4 in "ab":
d = ttx(r, o4)
dd = d[d.common_volts > 0]
rows.append(dict(o1=r[0], o2=r[1], o3=r[2], o4=o4, ttx_neurons=int((dd.group == "TTX").sum()), ttx_ok=int(dd[dd.group == "TTX"].ok.sum()),
vehicle_neurons=int((dd.group == "vehicle").sum()), vehicle_ok=int(dd[dd.group == "vehicle"].ok.sum()),
holds=bool(dd.ok.all()), ttx_ok_any=int(dd[dd.group == "TTX"].ok_any.sum()),
vehicle_ok_any=int(dd[dd.group == "vehicle"].ok_any.sum())))
TTXR = pd.DataFrame(rows); print(TTXR.to_string())
print(ttx(("c", "a", "a"), "a").to_string())
print("ok: every session at the common voltages responds before (TTX) or after (vehicle); ok_any: at least one does.")
print("longest pause (split) per vehicle neuron, s:", {n.split("/")[1]: round(g_) for n, g_ in ttx(("c", "a", "a"), "a").query("group == 'vehicle'")[["neuron", "gap_s"]].values})
o1 o2 o3 o4 ttx_neurons ttx_ok vehicle_neurons vehicle_ok holds ttx_ok_any vehicle_ok_any
0 c a a a 6 5 13 11 False 6 13
1 c a a b 6 5 13 11 False 6 13
2 c a b a 6 5 13 11 False 6 13
3 c a b b 6 5 13 11 False 6 13
4 c b a a 6 3 13 11 False 6 13
5 c b a b 6 3 13 11 False 6 13
6 c b b a 6 4 13 11 False 6 13
7 c b b b 6 4 13 11 False 6 13
8 a a a a 6 1 13 6 False 6 11
9 a a a b 6 1 13 6 False 6 11
10 a a b a 6 1 13 5 False 5 10
11 a a b b 6 1 13 5 False 5 10
12 a b a a 6 1 13 5 False 6 11
13 a b a b 6 1 13 5 False 6 11
14 a b b a 6 1 13 4 False 5 10
15 a b b b 6 1 13 4 False 5 10
16 b a a a 6 0 13 3 False 5 11
17 b a a b 6 0 13 3 False 5 11
18 b a b a 6 0 13 2 False 4 5
19 b a b b 6 0 13 2 False 4 5
20 b b a a 6 0 13 3 False 5 11
21 b b a b 6 0 13 3 False 5 11
22 b b b a 6 0 13 1 False 4 5
23 b b b b 6 0 13 1 False 4 5
neuron group gap_s common_volts before_resp after_resp aps_before aps_after sus_before sus_after ok ok_any
0 001271/10_20240731 TTX NaN 4 9/9 0/4 47.0 0.0 0.0 0.0 True True
1 001271/10_20240801 TTX NaN 6 7/7 0/6 80.0 0.0 0.0 0.0 True True
2 001271/12_20240731 vehicle 333.0 4 4/4 4/4 32.0 29.0 0.0 0.0 True True
3 001271/13_20240731 vehicle 167.0 4 5/5 4/4 37.0 27.0 0.0 0.0 True True
4 001271/13_20240801 vehicle 649.0 1 1/1 2/2 7.0 17.0 0.0 0.0 True True
5 001271/14_20240731 TTX NaN 5 7/8 0/5 62.0 0.0 2.0 0.0 False True
6 001271/14_20240801 vehicle 327.0 4 5/5 4/4 31.0 29.0 0.0 0.0 True True
7 001271/15_20240731 vehicle 312.0 5 5/5 5/5 54.0 50.0 0.0 0.0 True True
8 001271/15_20240801 vehicle 325.0 5 5/5 6/6 56.0 49.0 0.0 0.0 True True
9 001271/16_20240801 vehicle 372.0 7 6/7 7/7 50.0 40.0 0.0 0.0 True True
10 001271/17_20240801 vehicle 35.0 0 NaN NaN NaN NaN NaN NaN NaN NaN
11 001271/18_20240801 TTX NaN 6 6/6 0/6 37.0 0.0 0.0 0.0 True True
12 001271/2_20240730 TTX NaN 6 7/7 0/7 30.0 0.0 0.0 0.0 True True
13 001271/4_20240730 vehicle 73.0 3 4/4 2/3 7.0 5.0 0.0 0.0 False True
14 001271/5_20240730 TTX NaN 4 8/8 0/5 13.0 0.0 0.0 0.0 True True
15 001271/5_20240801 vehicle 321.0 7 6/7 5/7 37.0 27.0 0.0 0.0 False True
16 001271/6_20240801 vehicle 288.0 4 5/5 4/4 20.0 18.0 0.0 0.0 True True
17 001271/7_20240801 vehicle 282.0 4 4/4 4/4 5.0 10.0 0.0 0.0 True True
18 001271/8_20240801 vehicle 323.0 8 8/8 8/8 64.0 52.0 0.0 0.0 True True
19 001271/9_20240801 vehicle 301.0 6 7/7 6/6 53.0 45.0 0.0 0.0 True True
ok: every session at the common voltages responds before (TTX) or after (vehicle); ok_any: at least one does.
longest pause (split) per vehicle neuron, s: {'12_20240731': 333, '13_20240731': 167, '13_20240801': 649, '14_20240801': 327, '15_20240731': 312, '15_20240801': 325, '16_20240801': 372, '17_20240801': 35, '4_20240730': 73, '5_20240801': 321, '6_20240801': 288, '7_20240801': 282, '8_20240801': 323, '9_20240801': 301}
# The TTX sessions: mean dF in the last 100 ms of the stimulus against the O4 levels (ROI c, F0 a),
# and the neuron's sessions after TTX at its highest voltage.
g = rd(TT, ("c", "a", "a")); g = g[(g.kind == "sine") & g.ttx]
medh = rd(TT, ("c", "a", "a")).groupby("neuron").ap_h.apply(lambda s: np.median(np.concatenate(s.values)) if len(np.concatenate(s.values)) else np.nan)
s = g.groupby("neuron").agg(sessions=("path", "size"), vmax=("volts", "max"), last100_max=("mean_last100", "max"),
last100_median=("mean_last100", "median"), pre_sd=("pre_sd", "median"))
s["level_a"] = 5 * s.pre_sd; s["level_b"] = 0.5 * medh.reindex(s.index)
print(s.round(2).to_string())
# Per TTX neuron, at the voltages tested both before and after TTX: the mean dF in the last 100 ms of
# the stimulus (median over sessions), and the APs, before and after (ROI c, F0 a, O3 a).
G = rd(TT, ("c", "a", "a")); G = G[G.kind == "sine"]
rows = []
for n, x in G.groupby("neuron"):
if not x.ttx.any(): continue
b, a, _ = ts.split_before_after(x); common = sorted(set(b.volts) & set(a.volts))
for v in common:
bb, aa = b[b.volts == v], a[a.volts == v]
rows.append(dict(neuron=n, volts=v, last100_before=bb.mean_last100.median(), last100_after=aa.mean_last100.median(),
aps_before=int(bb.n_ap.sum()), aps_after=int(aa.n_ap.sum()), above5_before=int((bb.mean_last100 > 5).sum()), above5_after=int((aa.mean_last100 > 5).sum())))
print(pd.DataFrame(rows).round(2).to_string())
sessions vmax last100_max last100_median pre_sd level_a level_b
neuron
001271/10_20240731 10 12.0 0.42 -0.63 0.83 4.14 3.45
001271/10_20240801 26 6.0 -0.58 -0.69 0.79 3.97 4.25
001271/14_20240731 17 7.0 -0.19 -0.46 0.94 4.70 5.92
001271/18_20240801 6 6.5 -0.56 -0.74 0.79 3.96 3.90
001271/2_20240730 7 6.0 -0.35 -0.54 0.87 4.35 2.99
001271/5_20240730 5 6.5 -0.13 -0.42 0.93 4.67 2.63
neuron volts last100_before last100_after aps_before aps_after above5_before above5_after
0 001271/10_20240731 3.0 0.25 -0.68 7 0 0 0
1 001271/10_20240731 3.5 1.24 -0.56 15 0 0 0
2 001271/10_20240731 4.0 2.10 -0.51 15 0 0 0
3 001271/10_20240731 4.5 2.48 -0.82 10 0 0 0
4 001271/10_20240801 2.0 -0.61 -0.78 1 0 0 0
5 001271/10_20240801 2.5 0.06 -0.60 1 0 0 0
6 001271/10_20240801 3.0 1.13 -0.64 11 0 0 0
7 001271/10_20240801 3.5 2.52 -0.60 20 0 0 0
8 001271/10_20240801 4.0 3.02 -0.71 43 0 0 0
9 001271/10_20240801 4.5 2.70 -0.61 4 0 0 0
10 001271/14_20240731 2.5 0.09 -0.68 3 0 0 0
11 001271/14_20240731 3.0 1.04 -0.55 1 0 0 0
12 001271/14_20240731 4.0 4.21 -0.50 16 0 0 0
13 001271/14_20240731 4.5 5.08 -0.46 17 0 1 0
14 001271/14_20240731 5.0 5.18 -0.66 25 0 3 0
15 001271/18_20240801 4.0 0.42 -0.79 6 0 0 0
16 001271/18_20240801 4.5 0.59 -0.56 7 0 0 0
17 001271/18_20240801 5.0 1.49 -0.73 8 0 0 0
18 001271/18_20240801 5.5 1.57 -0.60 8 0 0 0
19 001271/18_20240801 6.0 1.70 -0.78 5 0 0 0
20 001271/18_20240801 6.5 1.73 -0.76 3 0 0 0
21 001271/2_20240730 2.5 -0.50 -0.69 3 0 0 0
22 001271/2_20240730 3.0 -0.58 -0.63 4 0 0 0
23 001271/2_20240730 3.5 -0.02 -0.40 5 0 0 0
24 001271/2_20240730 4.0 0.44 -0.63 7 0 0 0
25 001271/2_20240730 5.0 1.21 -0.39 6 0 0 0
26 001271/2_20240730 6.0 1.66 -0.44 5 0 0 0
27 001271/5_20240730 4.0 -0.56 -0.42 2 0 0 0
28 001271/5_20240730 5.0 -0.33 -0.58 3 0 0 0
29 001271/5_20240730 6.0 0.48 -0.28 4 0 0 0
30 001271/5_20240730 6.5 0.89 -0.36 4 0 0 0
Result: partially reproduced. Under the crop ROI all 6 TTX neurons had APs before TTX at the voltages tested in both conditions and none responded after it (no AP and no sustained depolarization); after TTX the mean dF in the last 100 ms is at most 0.42 %. In 5 of the 6 (3 or 4 with F0 b) every session before TTX responded; in 001271/14_20240731 one of 8 did not. Sustained depolarization (O4) before TTX occurred only in that neuron (2 sessions; its median last-100-ms dF was 5.08 and 5.18 % at 4.5 and 5 V before TTX, and -0.46 and -0.66 % after), so the abolition of the depolarization rests on one neuron. 11 of 13 vehicle neurons responded at every common voltage after the pause, and all 13 at some. The two that did not at every voltage are 001271/4_20240730 (2 of 3), whose longest pause is 73 s, so no change of solution can be placed, and 001271/5_20240801 (5 of 7). The vehicle split is the longest pause, about 300 s in most neurons but 167 s in 001271/13_20240731 and 73 s in 001271/4_20240730 (35 s in 001271/17_20240801, which has no voltage in common). Under the Otsu ROIs at most 1 of the 6 TTX neurons responds at every voltage before TTX.
fig3-not-heating¶
Excitation by pure sine waves was not caused by heating; the adiabatic temperature rise was 2.5 C after 250 ms at 63 V/cm in NaCl, at most 1.24 C for the same stimulus in the low-conductivity solution, where it evoked sustained firing, and 4.1 C at 81 V/cm in NaCl.
The rise is conductivity times E squared times duration over the volumetric heat capacity (water, 4.18 J per cm3 per K, the claims file's standard choice), with E the peak or the RMS field. These are calculations from the paper's own numbers; the Figure 3A neuron is not in the dandisets (no NaGlu 250-ms sessions).
def rise(E_vcm, sigma_mScm, t=0.25, rc=4.18e6):
return sigma_mScm / 10 * (E_vcm * 100) ** 2 * t / rc
for name, E, s, exp in [("63 V/cm NaCl", 63, 16.4, 2.5), ("63 V/cm NaGlu", 63, 8.4, 1.24), ("81 V/cm NaCl", 81, 16.4, 4.1)]:
print(f"{name}: peak {rise(E, s):.2f} C, RMS {rise(E, s) / 2:.2f} C, paper {exp}")
print("ratio 81/63 in NaCl: computed %.3f, paper %.3f (%.1f %% apart)" % ((81 / 63) ** 2, 4.1 / 2.5, 100 * abs((81 / 63) ** 2 / (4.1 / 2.5) - 1)))
print("ratio NaGlu/NaCl at 63: computed %.3f, paper %.3f (%.1f %% apart)" % (8.4 / 16.4, 1.24 / 2.5, 100 * abs((8.4 / 16.4) / (1.24 / 2.5) - 1)))
print("heat capacity that gives 2.5 C at 63 V/cm, peak field: %.2f J/cm3/K; RMS: %.2f" % (rise(63, 16.4) * 4.18 / 2.5, rise(63, 16.4) / 2 * 4.18 / 2.5))
63 V/cm NaCl: peak 3.89 C, RMS 1.95 C, paper 2.5 63 V/cm NaGlu: peak 1.99 C, RMS 1.00 C, paper 1.24 81 V/cm NaCl: peak 6.44 C, RMS 3.22 C, paper 4.1 ratio 81/63 in NaCl: computed 1.653, paper 1.640 (0.8 % apart) ratio NaGlu/NaCl at 63: computed 0.512, paper 0.496 (3.3 % apart) heat capacity that gives 2.5 C at 63 V/cm, peak field: 6.51 J/cm3/K; RMS: 3.25
Result: blocked, methods underspecified. With water's heat capacity, the peak field gives 3.89 C and the RMS field 1.95 C at 63 V/cm in NaCl, neither 2.5 C. The paper's three values scale as the formula says (81/63 squared is 1.653 against the paper's 1.640, 0.8 % apart; 8.4/16.4 is 0.512 against 0.496, 3.3 % apart), so its numbers are consistent with each other but come from a heat capacity or field the paper does not give (2.5 C needs 6.51 J per cm3 per K with the peak field). The Figure 3A neuron is not in the dandisets.
Strobe imaging: files and experiments¶
001273 (2-us pulses from two orthogonal directions, Figure 4) and 001274 (pulses and single beats of modulated 2-kHz sine waves, Figures 5, 6 and S3) hold full 512 x 512 frames and no ROIs. Their metadata are read first; the frames are streamed file by file in each section and reduced to ROI traces as they load.
from concurrent.futures import ThreadPoolExecutor
import re
SM = {d: [r for r in ts.load_dandiset(d, with_data=False)] for d in ["001273", "001274"]} # metadata only
def strobe_info(r):
s = r["session_description"]
v = re.search(r"at ([\d,.]+)\s*V amplitude", s)
fm = re.search(r"modulated at ([\d.]+)\s*Hz", s)
pos = "up/bottom" if "above (Up)" in s else ("right/left" if "Right and Left" in s else "not specified")
kind = "pulse" if "rectangular pulse" in s else ("mod" if fm else "sine")
return dict(path=r["path"], subject=r["subject_id"], start=r["start"], n=r["shape"][0], kind=kind, fm=float(fm[1]) if fm else np.nan,
volts=float(v[1].replace(",", ".")) if v else np.nan, backward="backward manner" in s, position=pos,
ttx="presence of the sodium channel inhibitor" in s, asset_id=r["asset_id"], dandiset=r["dandiset"])
url_cache = {}
def url_of(d, aid):
if aid not in url_cache: url_cache[aid] = ts.s3_url(d, aid)
return url_cache[aid]
def read_prepped(row):
for k in range(4):
try:
x = ts.read_strobe(url_of(row.dandiset, row.asset_id)).astype(np.float32)
break
except Exception:
time.sleep(3 + 3 * k)
if row.backward: x = x[::-1]
med, good = ts.frame_scale(x)
return x, med, good
S3 = pd.DataFrame([strobe_info(r) for r in SM["001273"] if "error" not in r])
S4 = pd.DataFrame([strobe_info(r) for r in SM["001274"] if "error" not in r])
print("001273:", len(S3), "files; failed:", sum("error" in r for r in SM["001273"]))
print(S3.groupby(["subject", "position", "volts"]).agg(files=("path", "size"), backward=("backward", "sum"), frames=("n", "first")).to_string())
print("\n001274:", len(S4), "files; failed:", sum("error" in r for r in SM["001274"]))
print(S4.groupby(["subject", "kind", "fm", "ttx"], dropna=False).agg(files=("path", "size"), frames=("n", "first"), volts=("volts", "first")).to_string())
notes = {m for r in SM["001273"] if "error" not in r for m in re.findall(r"[^.]*actual measured[^.]*\.", r["session_description"])}
print("001273 notes on the delivered amplitude:", notes)
001273: 181 files; failed: 0
files backward frames
subject position volts
1_20240117 not specified 40.0 7 4 160
80.0 8 4 160
1_20240119 right/left 40.0 8 4 160
up/bottom 40.0 8 4 160
1_20240311 up/bottom 40.0 8 4 160
2_20240117 not specified 40.0 8 4 160
80.0 8 4 160
2_20240119 right/left 40.0 8 4 160
up/bottom 40.0 8 4 160
100.0 1 0 160
2_20240311 up/bottom 40.0 8 4 160
3_20240117 not specified 40.0 8 4 160
80.0 8 4 160
3_20240119 right/left 40.0 8 4 160
100.0 2 0 160
up/bottom 40.0 8 4 160
100.0 1 0 160
3_20240311 up/bottom 40.0 8 4 160
4_20240117 not specified 40.0 8 4 160
80.0 8 4 160
4_20240119 right/left 40.0 8 4 160
100.0 1 0 160
up/bottom 40.0 8 4 160
100.0 1 0 160
4_20240311 up/bottom 60.0 8 4 160
6_20240119 up/bottom 100.0 8 4 160
130.0 8 4 160
001274: 232 files; failed: 0
files frames volts
subject kind fm ttx
1_20240522 pulse NaN False 2 140 10.0
sine NaN False 1 300 10.0
1_20240523 pulse NaN False 1 140 10.0
1_20240824 mod 154.0 False 8 325 7.0
True 8 325 7.0
pulse NaN False 8 160 10.0
True 8 160 10.0
2_20240522 mod 100.0 False 5 500 5.0
pulse NaN False 2 140 10.0
sine NaN False 2 300 5.0
2_20240523 mod 100.0 False 4 500 7.0
154.0 False 4 325 7.0
pulse NaN False 1 140 10.0
3_20240522 mod 100.0 False 4 500 7.0
pulse NaN False 1 140 10.0
sine NaN False 1 300 7.0
3_20240523 mod 100.0 False 4 500 7.0
154.0 False 4 325 7.0
pulse NaN False 1 140 10.0
3_20240824 mod 154.0 True 8 325 7.0
pulse NaN True 8 160 10.0
4_20240522 mod 100.0 False 4 500 7.0
pulse NaN False 1 140 15.0
sine NaN False 3 300 7.0
4_20240523 mod 100.0 False 4 500 7.0
154.0 False 3 325 7.0
pulse NaN False 1 140 10.0
4_20240824 mod 154.0 True 8 325 7.0
pulse NaN True 8 160 10.0
5_20240523 mod 100.0 False 4 500 7.0
154.0 False 4 325 7.0
pulse NaN False 1 140 10.0
5_20240824 mod 154.0 True 8 325 7.0
pulse NaN True 8 160 10.0
6_20240523 mod 100.0 False 4 500 7.0
154.0 False 4 325 7.0
pulse NaN False 1 140 10.0
6_20240824 mod 154.0 False 8 325 7.0
True 8 325 7.0
pulse NaN True 8 160 10.0
7_20240523 mod 100.0 False 4 500 7.0
154.0 False 4 325 7.0
pulse NaN False 1 140 10.0
7_20240824 mod 154.0 True 8 325 7.0
pulse NaN True 8 160 10.0
8_20240523 mod 100.0 False 4 500 7.0
154.0 False 4 325 7.0
8_20240824 mod 154.0 True 8 325 7.0
pulse NaN True 8 160 10.0
9_20240523 mod 100.0 False 4 500 7.0
154.0 False 4 325 7.0
001273 notes on the delivered amplitude: {' Note: For 40V pulse amplitude, the actual measured amplitude was 25V; for 80V, the actual measured was 62V; for 100V, the actual measured was 79V; and for 130V, the actual measured was 102V.'}
The paper's 19 Figure 4 experiments are cells in one field direction. 001273 has 22 groups of 7 or 8 replicates: on 2024-01-17 four cells at 40 V and again at 80 V with "electrode position not specified", on 2024-01-19 four cells with the electrodes up/bottom and right/left at 40 V and a fifth at 100 and 130 V, and on 2024-03-11 four cells up/bottom (one at 60 V). The paper gives 140 V/cm for all; which 19 groups are the paper's, and which cells are elongated or random, the files do not say. Two selections are run: (i) all 22 groups, (ii) the 15 groups at 40 V. In 001274 nine cells have a pulse and 2-kHz waves modulated at 100 Hz without TTX (2024-05-22 cells 2 to 4 and 2024-05-23 cells 2 to 7): these are taken as the nine cells of Figure 6. The paper does not mention TTX for these experiments, but 2024-08-24 sessions are labelled with TTX. The paper averaged 12 sine sessions for the Figure 5 cell; no cell in 001274 has more than 5 sessions at 100 Hz (2_20240522 has 5), so the Figure 5 and S3 cells (C7 and C8) cannot be identified.
fig4-pooled-time-constants¶
Membrane charging by 2-us pulses is fast; pooled over 19 experiments, charging and discharging time constants averaged 0.7 and 0.91 us at the cathodic pole and 0.65 and 0.8 us at the anodic pole, so the membrane charges fully within each sine wave cycle up to 100 to 200 kHz.
Per experiment, the replicates are reversed where recorded backward, scaled (reading a) or not (b), and averaged; ROIs from the averaged scaled stack. The cathodic pole is the depolarized ROI, the anodic the hyperpolarized one. Charging fit over the pulse (frames 20 to 60), discharging fit from the pulse end to the last frame with (a) a free amplitude, (b) the amplitude fixed to the charging fit's value at the pulse end. Expected: 0.7 +- 0.08, 0.91 +- 0.08, 0.65 +- 0.05, 0.8 +- 0.07 us.
groups = S3.groupby(["subject", "position", "volts"])
groups = [(k, g) for k, g in groups if len(g) >= 7]
print(len(groups), "experiments of 7 or more replicates")
tt = (np.arange(160) - ts.PULSE_ON) * 0.05
rows = []; t0 = time.time()
for (subj, pos, volts), g in groups:
with ThreadPoolExecutor(8) as ex:
reps = list(ex.map(read_prepped, [r for r in g.itertuples()]))
acc = {"a": np.zeros((160, 512, 512), np.float32), "b": np.zeros((160, 512, 512), np.float32)}; cnt = np.zeros(160)
for x, med, good in reps:
acc["a"][good] += x[good] / med[good][:, None, None]; acc["b"][good] += x[good]; cnt += good
goodall = cnt > 0
for k in acc: acc[k][goodall] /= cnt[goodall][:, None, None]
c1, c2, ratio, cell = ts.find_rois(acc["a"], goodall)
for o5 in "ab":
for pole, c in [("anodic", c1), ("cathodic", c2)]:
dF = ts.pulse_dF(ts.box(acc[o5], c), goodall, slice(1, ts.PULSE_ON))
A, tau_c = ts.fit_charging(tt[20:61], dF[20:61])
for rdg in "ab":
Afix = None if rdg == "a" else A * (1 - np.exp(-2.0 / tau_c))
A2, tau_d = ts.fit_discharging(tt[60:] - 2.0, dF[60:], Afix)
rows.append(dict(subject=subj, position=pos, volts=volts, replicates=len(g), o5=o5, pole=pole, fit=rdg,
dFmax=A, tau_charge=tau_c, tau_discharge=tau_d, ratio_at_roi=ratio[c]))
del acc, reps
print("read and fitted in %.0f s" % (time.time() - t0))
K4 = pd.DataFrame(rows)
print(K4[(K4.o5 == "a") & (K4.fit == "a")].round(2).to_string())
22 experiments of 7 or more replicates
read and fitted in 466 s
subject position volts replicates o5 pole fit dFmax tau_charge tau_discharge ratio_at_roi
0 1_20240117 not specified 40.0 7 a anodic a -7.64 0.71 1.59 0.92
2 1_20240117 not specified 40.0 7 a cathodic a 7.72 0.59 0.62 1.06
8 1_20240117 not specified 80.0 8 a anodic a -11.00 0.49 1.88 0.91
10 1_20240117 not specified 80.0 8 a cathodic a 38.60 8.10 2.00 1.08
16 1_20240119 right/left 40.0 8 a anodic a -4.34 0.31 2.72 0.97
18 1_20240119 right/left 40.0 8 a cathodic a 11.39 4.38 0.51 1.05
24 1_20240119 up/bottom 40.0 8 a anodic a -10.51 0.51 2.27 0.91
26 1_20240119 up/bottom 40.0 8 a cathodic a 20.27 1.13 1.90 1.14
32 1_20240311 up/bottom 40.0 8 a anodic a -16.44 0.80 0.67 0.86
34 1_20240311 up/bottom 40.0 8 a cathodic a 9.62 0.31 1.16 1.09
40 2_20240117 not specified 40.0 8 a anodic a -8.26 0.32 0.46 0.92
42 2_20240117 not specified 40.0 8 a cathodic a 8.86 0.74 0.44 1.08
48 2_20240117 not specified 80.0 8 a anodic a -13.83 0.46 0.47 0.87
50 2_20240117 not specified 80.0 8 a cathodic a 15.95 1.01 3.41 1.11
56 2_20240119 right/left 40.0 8 a anodic a -6.82 0.58 4.27 0.95
58 2_20240119 right/left 40.0 8 a cathodic a 10.51 2.11 1.00 1.07
64 2_20240119 up/bottom 40.0 8 a anodic a -6.51 1.68 0.43 0.96
66 2_20240119 up/bottom 40.0 8 a cathodic a 6.39 0.43 2.66 1.06
72 2_20240311 up/bottom 40.0 8 a anodic a -9.42 1.11 0.92 0.93
74 2_20240311 up/bottom 40.0 8 a cathodic a 5.16 0.38 2.52 1.05
80 3_20240117 not specified 40.0 8 a anodic a -9.43 0.49 0.66 0.90
82 3_20240117 not specified 40.0 8 a cathodic a 7.94 0.59 2.53 1.07
88 3_20240117 not specified 80.0 8 a anodic a -17.07 0.57 1.09 0.84
90 3_20240117 not specified 80.0 8 a cathodic a 7.59 0.45 1.08 1.07
96 3_20240119 right/left 40.0 8 a anodic a -64.56 14.36 3.39 0.94
98 3_20240119 right/left 40.0 8 a cathodic a 4.46 0.19 0.01 1.06
104 3_20240119 up/bottom 40.0 8 a anodic a -81.83 20.00 1.59 0.94
106 3_20240119 up/bottom 40.0 8 a cathodic a 7.94 0.31 1.07 1.08
112 3_20240311 up/bottom 40.0 8 a anodic a -98.24 20.00 2.17 0.91
114 3_20240311 up/bottom 40.0 8 a cathodic a 7.38 0.22 1.99 1.06
120 4_20240117 not specified 40.0 8 a anodic a -10.31 1.03 0.35 0.91
122 4_20240117 not specified 40.0 8 a cathodic a 8.44 0.49 2.20 1.07
128 4_20240117 not specified 80.0 8 a anodic a -11.38 0.57 1.65 0.90
130 4_20240117 not specified 80.0 8 a cathodic a 7.49 0.43 0.90 1.07
136 4_20240119 right/left 40.0 8 a anodic a -5.88 0.36 0.20 0.93
138 4_20240119 right/left 40.0 8 a cathodic a 9.62 0.54 1.59 1.09
144 4_20240119 up/bottom 40.0 8 a anodic a -95.17 20.00 0.38 0.92
146 4_20240119 up/bottom 40.0 8 a cathodic a 5.66 0.16 1.70 1.05
152 4_20240311 up/bottom 60.0 8 a anodic a -13.44 0.54 1.93 0.88
154 4_20240311 up/bottom 60.0 8 a cathodic a 14.65 0.80 0.88 1.13
160 6_20240119 up/bottom 100.0 8 a anodic a -28.76 0.69 1.26 0.74
162 6_20240119 up/bottom 100.0 8 a cathodic a 14.99 0.40 1.32 1.14
168 6_20240119 up/bottom 130.0 8 a anodic a -24.62 0.58 1.62 0.78
170 6_20240119 up/bottom 130.0 8 a cathodic a 11.51 0.44 1.29 1.11
def pooled(df):
out = {}
for pole in ["cathodic", "anodic"]:
x = df[df.pole == pole]
out[f"{pole}_charge"] = f"{x.tau_charge.mean():.2f} +- {x.tau_charge.sem():.2f}"
out[f"{pole}_discharge"] = f"{x.tau_discharge.mean():.2f} +- {x.tau_discharge.sem():.2f}"
out[f"{pole}_charge_median"] = x.tau_charge.median(); out[f"{pole}_discharge_median"] = x.tau_discharge.median()
exp = {"cathodic_charge": (0.7, 0.08), "cathodic_discharge": (0.91, 0.08), "anodic_charge": (0.65, 0.05), "anodic_discharge": (0.8, 0.07)}
out["within_sem"] = sum(abs(float(out[k].split()[0]) - v) <= tol for k, (v, tol) in exp.items())
return out
rows = []
for sel, f in [("all 22", lambda d: d), ("40 V only", lambda d: d[d.volts == 40])]:
for o5 in "ab":
for fit in "ab":
d = f(K4[(K4.o5 == o5) & (K4.fit == fit)])
rows.append(dict(selection=sel, o5=o5, fit=fit, experiments=len(d) // 2, **pooled(d)))
P4 = pd.DataFrame(rows)
print(P4.to_string())
print("time to 95%% charge (3 tau) for the largest charging tau over experiments: %.2f us" % (3 * K4.tau_charge.max()))
selection o5 fit experiments cathodic_charge cathodic_discharge cathodic_charge_median cathodic_discharge_median anodic_charge anodic_discharge anodic_charge_median anodic_discharge_median within_sem 0 all 22 a a 22 1.10 +- 0.39 1.49 +- 0.18 0.469866 1.301194 3.92 +- 1.53 1.45 +- 0.23 0.576650 1.423464 0 1 all 22 a b 22 1.10 +- 0.39 1.42 +- 0.20 0.469866 1.138567 3.92 +- 1.53 1.09 +- 0.13 0.576650 1.065720 0 2 all 22 b a 22 1.86 +- 0.94 3.18 +- 1.18 0.463658 1.436049 4.30 +- 1.63 2.53 +- 0.84 0.566174 1.741296 0 3 all 22 b b 22 1.86 +- 0.94 1.46 +- 0.25 0.463658 1.235760 4.30 +- 1.63 1.39 +- 0.22 0.566174 1.130930 0 4 40 V only a a 15 0.84 +- 0.28 1.46 +- 0.22 0.486195 1.587116 5.48 +- 2.14 1.47 +- 0.32 0.797498 0.916431 0 5 40 V only a b 15 0.84 +- 0.28 1.32 +- 0.25 0.486195 1.102115 5.48 +- 2.14 1.05 +- 0.18 0.797498 0.813198 0 6 40 V only b a 15 1.18 +- 0.54 4.02 +- 1.69 0.368464 1.786939 6.04 +- 2.27 2.93 +- 1.23 1.024211 1.696103 0 7 40 V only b b 15 1.18 +- 0.54 1.40 +- 0.30 0.368464 1.351222 6.04 +- 2.27 1.41 +- 0.31 1.024211 0.989254 0 time to 95% charge (3 tau) for the largest charging tau over experiments: 60.00 us
# Not deciding: pole and kinetic comparisons over all 22 experiments (scaled, fit a), for the blocked
# fig4-orientation-and-poles. Paired two-sided t tests.
x = K4[(K4.o5 == "a") & (K4.fit == "a")].pivot_table(index=["subject", "position", "volts"], columns="pole", values=["dFmax", "tau_charge", "tau_discharge"])
for v in ["dFmax", "tau_charge", "tau_discharge"]:
a, c = x[(v, "anodic")].abs(), x[(v, "cathodic")].abs()
print(f"{v}: anodic median {a.median():.2f}, cathodic median {c.median():.2f}, paired p = {stats.ttest_rel(a, c).pvalue:.3g}")
for pole in ["anodic", "cathodic"]:
print(f"{pole}: discharging minus charging tau median {(x[('tau_discharge', pole)] - x[('tau_charge', pole)]).median():.2f} us, paired p = {stats.ttest_rel(x[('tau_discharge', pole)], x[('tau_charge', pole)]).pvalue:.3g}")
print("charging fits at the 20-us bound:", int((K4.tau_charge >= 19.99).sum()), "of", len(K4))
dFmax: anodic median 11.19, cathodic median 8.65, paired p = 0.0533 tau_charge: anodic median 0.58, cathodic median 0.47, paired p = 0.105 tau_discharge: anodic median 1.42, cathodic median 1.30, paired p = 0.917 anodic: discharging minus charging tau median 0.15 us, paired p = 0.119 cathodic: discharging minus charging tau median 0.80 us, paired p = 0.368 charging fits at the 20-us bound: 16 of 176
Result: blocked, metadata insufficient. The claim is a mean over the paper's 19 experiments, measured in hand-drawn ROIs. 001273 holds no ROIs and does not say which of its 22 groups are the 19 experiments; its groups were run at 40 to 130 V, whose session descriptions give measured amplitudes of 25 to 102 V, where the paper gives 140 V/cm for all. Shown, not deciding: with 12 x 12 pixel boxes at the extremes of the smoothed ratio image, no pooled mean is within the paper's SEM under any selection, scaling or discharging fit (0 of 4 in every row), 16 of 176 charging fits end at the 20-us bound, and the medians are 0.37 to 0.49 us (cathodic) and 0.57 to 1.02 us (anodic) for charging and 0.81 to 1.79 us for discharging, against 0.65 to 0.91 in the paper. The first round recorded this claim as not reproduced; the time constants depend on the ROIs and the selection, which the dandiset does not give.
fig6-sine-less-than-pulse¶
Across 9 cells FluoVolt sensitivity to 2-us pulses was 15% to 47% per 1 kV/cm (mean 28%); sensitivity to 2-kHz sine waves, linearly correlated with the field, was 1.5% to 28.4% (mean 11.8%), smaller in every ROI, and 42% of the pulse value.
Per cell: ROIs from the pulse (above); pulse sensitivity, the mean dF over the 21 frames from 1 to 2 us into the pulse over 0.257 kV/cm, sign removed. Sine: the 2-kHz waves modulated at 100 Hz (four sessions per cell, five in 2_20240522), dF per session (F0 from the first and last cycles), averaged over sessions, the LOESS fit subtracted (R's loess, degree 2, as in RStudio, Table S1; O6: span 2 or 4 ms), and the oscillation regressed (a) on the signed field (sign removed after), (b) per carrier cycle, half the peak-to-peak amplitude on the cycle's peak field (half, so that it is on the scale of (a); the criterion names the full peak-to-peak amplitude, shown as a variant). Scaling O5 (a) or (b). As a check of the assumed carrier phase, the phase of the carrier in each ROI's oscillation is also fitted, with the sensitivity at that phase (variant c).
cells6 = []
for subj, g in S4[~S4.ttx].groupby("subject"):
if (g.kind == "pulse").any() and ((g.kind == "mod") & (g.fm == 100)).any():
cells6.append(subj)
print(len(cells6), "cells:", cells6)
def read_rois(row, centres):
for k in range(4):
try:
x = ts.read_strobe(url_of(row.dandiset, row.asset_id)); break
except Exception:
time.sleep(3 + 3 * k)
if row.backward: x = x[::-1]
med, good = ts.frame_scale(x)
tr = np.array([[ts.box(x, c, h) for c in centres] for h in ROI_HALVES])
return tr, med, good
ROI_HALVES = [6, 3, 10] # 12 x 12 boxes used; 6 x 6 and 20 x 20 as a variant (round-2 verifier: the sizes depend on the ROIs)
SENS = []; WHOLE = []; PULSE_TR = []; BOXV = []; t0 = time.time()
for subj in cells6:
g = S4[(S4.subject == subj) & ~S4.ttx]
pr = g[g.kind == "pulse"].sort_values("start").iloc[0]
x, med, good = read_prepped(pr)
xs = x / med[:, None, None]
c1, c2, ratio, cell = ts.find_rois(xs, good)
sess = list(g[(g.kind == "mod") & (g.fm == 100)].itertuples())
with ThreadPoolExecutor(4) as ex:
sine = list(ex.map(lambda r: read_rois(r, (c1, c2)), sess))
n = min(s[0].shape[-1] for s in sine)
v_sine = g[(g.kind == "mod") & (g.fm == 100)].volts.iloc[0]
epk = 0.136 * v_sine / 7.0 # 136 V/cm at the 7 V of the Figure 5 cell; scaled for cells run at another voltage
epulse = 0.257 * pr.volts / 10.0 # 257 V/cm at 10 V, likewise
t, E, env = ts.sine_field(n, 100.0, epk=epk)
base = np.r_[0:25, n - 25:n]
for o5 in "ab":
P = xs if o5 == "a" else x
for k, (roi, c) in enumerate([("ROI1", c1), ("ROI2", c2)]):
dFp = ts.pulse_dF(ts.box(P, c), good, slice(1, ts.PULSE_ON))
ps = np.nanmean(dFp[40:61]) / epulse
run5 = pd.Series(dFp[ts.PULSE_ON:ts.PULSE_OFF + 1]).rolling(5, center=True).mean()
peak5 = run5.abs().max() / epulse
dfs = []
for tr, m, gd in sine:
F = tr[0, k, :n] / (m[:n] if o5 == "a" else 1.0); F = F.astype(float); F[~gd[:n]] = np.nan
dfs.append(100 * (F / np.nanmean(F[base]) - 1))
dF = np.nanmean(dfs, 0); ok = np.isfinite(dF)
for span in [2.0, 4.0]:
lo = ts.loess(t[ok], dF[ok], span / (n * 0.02))
osc = dF[ok] - lo
lr = stats.linregress(E[ok], osc)
cyc = (t[ok] // 5e-4).astype(int)
cp = pd.DataFrame(dict(c=cyc, o=osc, e=env[ok] * epk)).groupby("c").agg(p2p=("o", lambda v: v.max() - v.min()), epk=("e", "max"))
lr_b = stats.linregress(cp.epk, cp.p2p / 2)
phase, amp = ts.carrier_phase(t[ok], osc, env[ok])
SENS.append(dict(cell=subj, roi=roi, o5=o5, span=span, pulse=abs(ps), pulse_peak5=peak5, pulse_signed=ps,
sine_a=abs(lr.slope), sine_a_signed=lr.slope, p_a=lr.pvalue, sine_b=abs(lr_b.slope), p_b=lr_b.pvalue,
sessions=len(sine), sine_volts=v_sine, pulse_volts=pr.volts, phase=phase, sine_c=amp / epk, pulse_frames=pr.n))
if o5 == "a" and span == 2.0:
PULSE_TR.append(dFp * np.sign(ps))
if o5 == "a":
WHOLE.append(dict(cell=subj, roi=roi, span=span, loess_at_peak=lo[np.argmax(env[ok])], osc=osc))
# Box-size variant: scaled (O5 a), span 2 ms, signed regression, same centres.
for hi, h in enumerate(ROI_HALVES[1:], 1):
for k, c in enumerate((c1, c2)):
dFp = ts.pulse_dF(ts.box(xs, c, h), good, slice(1, ts.PULSE_ON))
ps = np.nanmean(dFp[40:61]) / epulse
dfs = []
for tr, m, gd in sine:
F = (tr[hi, k, :n] / m[:n]).astype(float); F[~gd[:n]] = np.nan
dfs.append(100 * (F / np.nanmean(F[base]) - 1))
dF = np.nanmean(dfs, 0); ok = np.isfinite(dF)
osc = dF[ok] - ts.loess(t[ok], dF[ok], 2.0 / (n * 0.02))
lr = stats.linregress(E[ok], osc)
BOXV.append(dict(cell=subj, roi=f"ROI{k + 1}", half=h, pulse=abs(ps), sine_a=abs(lr.slope), p_a=lr.pvalue))
print("read in %.0f s" % (time.time() - t0))
SE = pd.DataFrame(SENS)
print(SE[(SE.o5 == "a") & (SE.span == 2.0)].drop(columns=["o5", "span"]).round(4).to_string())
# The pulse onset: the files say exposures began 1 us after the start of recording for single pulses
# (frame 20 at 50-ns steps), for the 140- and the 160-frame files alike. The ROI dF (sign set so that
# the pulse response is positive, scaled, mean over the 18 ROIs) in 0.25-us bins from the first frame:
print("frames of the pulse files used:", sorted(SE.pulse_frames.unique().tolist()))
pt_ = np.nanmean(np.array([p[:140] for p in PULSE_TR]), 0)
print("mean dF per 0.25-us bin (5 frames) from frame 0:", np.round([np.nanmean(pt_[k:k + 5]) for k in range(0, 140, 5)], 1).tolist())
9 cells: ['2_20240522', '2_20240523', '3_20240522', '3_20240523', '4_20240522', '4_20240523', '5_20240523', '6_20240523', '7_20240523']
read in 239 s
cell roi pulse pulse_peak5 pulse_signed sine_a sine_a_signed p_a sine_b p_b sessions sine_volts pulse_volts phase sine_c pulse_frames
0 2_20240522 ROI1 14.8352 20.4376 -14.8352 0.7120 -0.7120 0.5308 5.4962 0.2340 5 5.0 10.0 -1.8202 2.8842 140
2 2_20240522 ROI2 16.6237 21.8014 16.6237 8.7899 8.7899 0.0000 4.0992 0.6317 5 5.0 10.0 0.4455 9.7405 140
8 2_20240523 ROI1 27.0897 33.4831 -27.0897 11.2640 -11.2640 0.0000 21.0226 0.2551 4 7.0 10.0 -3.1410 11.2640 140
10 2_20240523 ROI2 8.3306 14.4781 8.3306 0.0350 0.0350 0.9775 19.4992 0.1613 4 7.0 10.0 1.5567 2.4849 140
16 3_20240522 ROI1 21.5308 25.1826 -21.5308 5.9587 -5.9587 0.0000 11.5161 0.3155 4 7.0 10.0 -2.6386 6.8012 140
18 3_20240522 ROI2 6.8767 13.1662 6.8767 3.2638 3.2638 0.0442 23.9675 0.2070 4 7.0 10.0 0.2374 3.3580 140
24 3_20240523 ROI1 16.4327 21.5736 -16.4327 6.7220 -6.7220 0.0000 1.6428 0.5336 4 7.0 10.0 -2.6317 7.7015 140
26 3_20240523 ROI2 13.8988 24.0369 13.8988 4.7125 4.7125 0.0000 0.6912 0.7178 4 7.0 10.0 0.1268 4.7506 140
32 4_20240522 ROI1 20.8021 25.3315 -20.8021 7.9649 -7.9649 0.0000 6.9936 0.0503 4 7.0 15.0 -2.7068 8.7822 140
34 4_20240522 ROI2 15.9303 19.2699 15.9303 6.8499 6.8499 0.0000 0.1387 0.9526 4 7.0 15.0 0.2617 7.0914 140
40 4_20240523 ROI1 7.4816 16.0070 -7.4816 0.9788 -0.9788 0.2812 4.0066 0.1230 4 7.0 10.0 -2.6957 1.0849 140
42 4_20240523 ROI2 3.1372 12.0666 3.1372 1.4280 -1.4280 0.1804 5.2886 0.2304 4 7.0 10.0 -2.1264 2.7075 140
48 5_20240523 ROI1 18.9701 22.5880 -18.9701 3.9313 -3.9313 0.0000 2.3405 0.3034 4 7.0 10.0 -2.3850 5.4060 140
50 5_20240523 ROI2 24.1916 31.2744 24.1916 4.1221 4.1221 0.0000 1.3908 0.3788 4 7.0 10.0 0.1109 4.1476 140
56 6_20240523 ROI1 29.1709 33.6261 -29.1709 4.0784 -4.0784 0.0000 2.1544 0.3027 4 7.0 10.0 -2.2228 6.7217 140
58 6_20240523 ROI2 23.3567 30.9956 23.3567 12.2746 12.2746 0.0000 4.8426 0.0058 4 7.0 10.0 0.1824 12.4817 140
64 7_20240523 ROI1 26.0151 30.7161 -26.0151 14.8096 -14.8096 0.0000 6.2513 0.0006 4 7.0 10.0 -2.9319 15.1411 140
66 7_20240523 ROI2 17.6877 26.8079 17.6877 2.0219 2.0219 0.0074 2.3303 0.1865 4 7.0 10.0 0.8860 3.1967 140
frames of the pulse files used: [140]
mean dF per 0.25-us bin (5 frames) from frame 0: [0.3, 0.3, -0.2, -0.0, 1.0, 2.9, 3.8, 4.2, 4.2, 4.3, 5.3, 5.3, 2.8, 1.3, 0.9, 0.4, 0.6, -0.0, -0.5, 0.6, -0.2, 0.0, -1.1, -1.0, -1.2, -1.2, -1.4, -1.5]
rows = []
for o5 in "ab":
for span in [2.0, 4.0]:
for sl in ["a", "b"]:
x = SE[(SE.o5 == o5) & (SE.span == span)]
s = x[f"sine_{sl}"]; p = x[f"p_{sl}"]
cells_sig = x.assign(sig=p < 1e-4).groupby("cell").sig.all().sum()
d = dict(o5=o5, span=span, slope=sl, rois=len(x), pulse_mean=x.pulse.mean(), pulse_sem=x.pulse.sem(), pulse_min=x.pulse.min(),
pulse_max=x.pulse.max(), sine_mean=s.mean(), sine_sem=s.sem(), sine_min=s.min(), sine_max=s.max(),
ratio_pct=100 * s.mean() / x.pulse.mean(), sine_below_pulse=int((s < x.pulse).sum()), cells_p_both=int(cells_sig))
d["holds"] = (abs(d["pulse_mean"] - 28) <= 2.5 and abs(d["sine_mean"] - 11.8) <= 1.7 and abs(d["ratio_pct"] - 42) <= 4.2
and 13.5 <= d["pulse_min"] and d["pulse_max"] <= 51.7 and 1.35 <= d["sine_min"] and d["sine_max"] <= 31.24
and d["sine_below_pulse"] == d["rois"] and abs(d["cells_p_both"] - 8) <= 0.8)
rows.append(d)
F6 = pd.DataFrame(rows)
with pd.option_context("display.float_format", "{:.3g}".format):
print(F6.to_string())
# Reading (b) with the full peak-to-peak amplitude, as the criterion names it (round-3 verifier: the
# table above regresses half of it); the slopes double and the p-values are the same.
for o5 in "ab":
for span in [2.0, 4.0]:
x = SE[(SE.o5 == o5) & (SE.span == span)]; s = 2 * x.sine_b
h = (abs(x.pulse.mean() - 28) <= 2.5 and abs(s.mean() - 11.8) <= 1.7 and abs(100 * s.mean() / x.pulse.mean() - 42) <= 4.2
and 13.5 <= x.pulse.min() and x.pulse.max() <= 51.7 and 1.35 <= s.min() and s.max() <= 31.24
and (s < x.pulse).all() and abs(x.assign(sig=x.p_b < 1e-4).groupby("cell").sig.all().sum() - 8) <= 0.8)
print(f"reading b, full peak-to-peak, o5 {o5}, span {span}: sine mean {s.mean():.2f}, range {s.min():.2f} to {s.max():.2f}, "
f"ratio {100 * s.mean() / x.pulse.mean():.1f} %, below pulse in {int((s < x.pulse).sum())} of {len(x)}, holds {h}")
x = SE[(SE.o5 == "a") & (SE.span == 2.0)]
print("antiphase: sign of the signed sine slope in ROI1 / ROI2 per cell:",
{c: tuple(np.sign(v.sine_a_signed.values).astype(int)) for c, v in x.groupby("cell")})
# The carrier phase (round-1 verifier, problem 10): the phase of the carrier in each ROI's oscillation,
# fitted freely (0 = in phase with the assumed field, +-pi = antiphase), and the sensitivity at the
# fitted phase (amplitude over peak field), which does not rest on the assumed phase (variant c).
print("fitted carrier phase (rad) ROI1 / ROI2 per cell:", {c: tuple(np.round(v.phase.values, 2)) for c, v in x.groupby("cell")})
for o5 in "ab":
for span in [2.0, 4.0]:
y = SE[(SE.o5 == o5) & (SE.span == span)]
print(f"variant c (fitted phase), o5 {o5}, span {span}: sine mean {y.sine_c.mean():.2f} +- {y.sine_c.sem():.2f}, range {y.sine_c.min():.2f} to {y.sine_c.max():.2f}, "
f"ratio {100 * y.sine_c.mean() / y.pulse.mean():.1f} %, below pulse in {int((y.sine_c < y.pulse).sum())} of {len(y)}")
print("pulse: sign of ROI1 / ROI2:", {c: tuple(np.sign(v.pulse_signed.values).astype(int)) for c, v in x.groupby("cell")})
o5 span slope rois pulse_mean pulse_sem pulse_min pulse_max sine_mean sine_sem sine_min sine_max ratio_pct sine_below_pulse cells_p_both holds
0 a 2 a 18 17.4 1.75 3.14 29.2 5.55 0.99 0.035 14.8 32 18 4 False
1 a 2 b 18 17.4 1.75 3.14 29.2 6.87 1.72 0.139 24 39.6 15 0 False
2 a 4 a 18 17.4 1.75 3.14 29.2 5.61 1 0.0213 15 32.3 18 4 False
3 a 4 b 18 17.4 1.75 3.14 29.2 7.03 1.81 0.0899 24.4 40.5 15 0 False
4 b 2 a 18 17.4 1.7 4.72 28.1 5.55 1.01 0.609 15.8 31.9 18 4 False
5 b 2 b 18 17.4 1.7 4.72 28.1 7.8 2.23 0.0441 27 44.8 15 0 False
6 b 4 a 18 17.4 1.7 4.72 28.1 5.61 1.02 0.598 16 32.2 18 4 False
7 b 4 b 18 17.4 1.7 4.72 28.1 8.08 2.35 0.0564 27.9 46.4 15 0 False
reading b, full peak-to-peak, o5 a, span 2.0: sine mean 13.74, range 0.28 to 47.94, ratio 79.2 %, below pulse in 12 of 18, holds False
reading b, full peak-to-peak, o5 a, span 4.0: sine mean 14.05, range 0.18 to 48.76, ratio 81.0 %, below pulse in 13 of 18, holds False
reading b, full peak-to-peak, o5 b, span 2.0: sine mean 15.61, range 0.09 to 54.02, ratio 89.6 %, below pulse in 12 of 18, holds False
reading b, full peak-to-peak, o5 b, span 4.0: sine mean 16.17, range 0.11 to 55.84, ratio 92.8 %, below pulse in 13 of 18, holds False
antiphase: sign of the signed sine slope in ROI1 / ROI2 per cell: {'2_20240522': (np.int64(-1), np.int64(1)), '2_20240523': (np.int64(-1), np.int64(1)), '3_20240522': (np.int64(-1), np.int64(1)), '3_20240523': (np.int64(-1), np.int64(1)), '4_20240522': (np.int64(-1), np.int64(1)), '4_20240523': (np.int64(-1), np.int64(-1)), '5_20240523': (np.int64(-1), np.int64(1)), '6_20240523': (np.int64(-1), np.int64(1)), '7_20240523': (np.int64(-1), np.int64(1))}
fitted carrier phase (rad) ROI1 / ROI2 per cell: {'2_20240522': (np.float64(-1.82), np.float64(0.45)), '2_20240523': (np.float64(-3.14), np.float64(1.56)), '3_20240522': (np.float64(-2.64), np.float64(0.24)), '3_20240523': (np.float64(-2.63), np.float64(0.13)), '4_20240522': (np.float64(-2.71), np.float64(0.26)), '4_20240523': (np.float64(-2.7), np.float64(-2.13)), '5_20240523': (np.float64(-2.39), np.float64(0.11)), '6_20240523': (np.float64(-2.22), np.float64(0.18)), '7_20240523': (np.float64(-2.93), np.float64(0.89))}
variant c (fitted phase), o5 a, span 2.0: sine mean 6.43 +- 0.91, range 1.08 to 15.14, ratio 37.1 %, below pulse in 18 of 18
variant c (fitted phase), o5 a, span 4.0: sine mean 6.50 +- 0.92, range 1.10 to 15.32, ratio 37.5 %, below pulse in 18 of 18
variant c (fitted phase), o5 b, span 2.0: sine mean 6.40 +- 0.94, range 1.00 to 16.14, ratio 36.7 %, below pulse in 18 of 18
variant c (fitted phase), o5 b, span 4.0: sine mean 6.47 +- 0.95, range 1.02 to 16.33, ratio 37.1 %, below pulse in 18 of 18
pulse: sign of ROI1 / ROI2: {'2_20240522': (np.int64(-1), np.int64(1)), '2_20240523': (np.int64(-1), np.int64(1)), '3_20240522': (np.int64(-1), np.int64(1)), '3_20240523': (np.int64(-1), np.int64(1)), '4_20240522': (np.int64(-1), np.int64(1)), '4_20240523': (np.int64(-1), np.int64(1)), '5_20240523': (np.int64(-1), np.int64(1)), '6_20240523': (np.int64(-1), np.int64(1)), '7_20240523': (np.int64(-1), np.int64(1))}
# Box-size variant (round-2 verifier): the 12 x 12 boxes are a choice the files do not settle, and
# the paper gives no ROI size. The same centres with 6 x 6 and 20 x 20 boxes (scaled, span 2 ms,
# signed regression).
B = pd.concat([pd.DataFrame(BOXV), SE[(SE.o5 == "a") & (SE.span == 2.0)][["cell", "roi", "pulse", "sine_a", "p_a"]].assign(half=6)], ignore_index=True)
rows = []
for h, x in B.groupby("half"):
rows.append(dict(box=f"{2 * h} x {2 * h}", rois=len(x), pulse_mean=x.pulse.mean(), pulse_sem=x.pulse.sem(), pulse_min=x.pulse.min(), pulse_max=x.pulse.max(),
sine_mean=x.sine_a.mean(), sine_min=x.sine_a.min(), sine_max=x.sine_a.max(), ratio_pct=100 * x.sine_a.mean() / x.pulse.mean(),
sine_below_pulse=int((x.sine_a < x.pulse).sum()), cells_p_both=int(x.assign(sig=x.p_a < 1e-4).groupby("cell").sig.all().sum())))
with pd.option_context("display.float_format", "{:.3g}".format):
print(pd.DataFrame(rows).to_string())
box rois pulse_mean pulse_sem pulse_min pulse_max sine_mean sine_min sine_max ratio_pct sine_below_pulse cells_p_both 0 6 x 6 18 16.9 1.96 4.51 30.2 5.96 0.607 14 35.2 17 1 1 12 x 12 18 17.4 1.75 3.14 29.2 5.55 0.035 14.8 32 18 4 2 20 x 20 18 13.3 1.45 3.66 22.4 4.98 0.00803 13.7 37.6 18 4
Result: blocked, metadata insufficient. Six of the criterion's seven parts are values (means, ranges, ratio, the count of cells with p < 0.0001) measured in the paper's hand-drawn ROIs on the paper's nine cells, and 001274 holds neither the ROIs nor labels for the cells. Both sensitivities are also divided by the paper's nominal fields (257 and 136 V/cm), scaled by the generator voltage, as the files give no measured amplitudes; a wrong field scales every value. Not deciding: sine sensitivity is below pulse sensitivity in all 18 ROIs of the 9 cells under the signed regression (a), and in 15 of 18 under the per-cycle reading (b). With the 12 x 12 boxes the pulse sensitivity averages 17.4 % per kV/cm (range 3.1 to 29.2; paper 28, 15 to 47), the sine sensitivity 5.55 % (0.04 to 14.8; paper 11.8, 1.5 to 28.4) and their ratio 32 % (paper 42 %), and the regression reaches p < 0.0001 in both ROIs in 4 cells, not 8. These values depend on the ROIs, which the paper drew by hand and the files do not hold: with 6 x 6 and 20 x 20 boxes at the same centres the pulse mean is 16.9 and 13.3 %, the sine mean 5.96 and 4.98 %, the ratio 35.2 and 37.6 %, and 1 and 4 cells reach p < 0.0001 in both ROIs; the sine value is below the pulse in 17 and 18 of 18 ROIs. The nine cells are also inferred, not labelled. As for fig4-pooled-time-constants, which was blocked because its values depend on ROIs and a selection the dandiset does not give, the sizes cannot be judged here; they are outside the paper's under every box size and reading, but that is not shown to be a failure of the paper. The carrier phase is assumed zero at the first frame. Fitted freely, the ROI 1 carriers lie at -1.8 to -3.1 rad in all 9 cells and the ROI 2 carriers at 0.1 to 0.9 rad in 7; the two ROIs are near antiphase in those 7, while 4_20240523 is nearly in phase (-2.70 and -2.13 rad) and 2_20240523 in quadrature (-3.14 and 1.56 rad). At the fitted phase the sine sensitivity averages 6.4 to 6.5 % (ratio 37 %), below the pulse value in 18 of 18 ROIs. The pulse files' description puts the pulse onset at frame 20 (1 us), and the mean pulse dF is near 0 before it, rises from frames 20 to 29 and falls after frame 60. The direction holds in every ROI under the signed regression and at the fitted phase, in 15 of 18 under the per-cycle reading, and in 17 or 18 of 18 under every box size; it is the only part that does not rest on the paper's ROIs, and as it fails under the per-cycle reading it cannot carry the claim. Under the per-cycle reading no cell reaches p < 0.0001 in both ROIs. With the full peak-to-peak amplitude, as the criterion names it, the per-cycle sine mean is 13.7 to 16.2 %, the ratio 79 to 93 %, and the sine value is below the pulse in 12 or 13 of 18 ROIs; no row holds.
threshold-depolarization-estimate¶
The 2-kHz AP threshold of 57.7 V/cm times the sine wave sensitivity gives a 0.68% emission change, about 6.8 mV with FluoVolt's calibration, so the threshold depolarization is 5 to 15 mV, as with intracellular current injection.
The mean sine sensitivity of fig6-sine-less-than-pulse times the mean unmodulated 2-kHz threshold of fig2d-modulated-thresholds-higher (1-AP or 5-AP, under both field readings), over 1,000, and that times 10 mV per 1%.
rows = []
for r in [("c", "a", "a"), ("a", "a", "a"), ("b", "a", "a")]:
th = TH[r]; x = th[th.stim == "pure_2k"]
for k in ["t1", "t5"]:
for fr, fac in ts.FIELD_READINGS.items():
thr = x[k].mean() * fac * ts.V_TO_VCM
for o5 in "ab":
for span in [2.0, 4.0]:
for sl in "ab":
s = SE[(SE.o5 == o5) & (SE.span == span)][f"sine_{sl}"].mean()
em = s * thr / 1000
rows.append(dict(roi_reading="".join(r), threshold=k, field=fr, threshold_vcm=thr, o5=o5, span=span, slope=sl,
sine=s, emission_pct=em, mV=em * 10, ok=abs(em - 0.68) <= 0.11 and abs(em * 10 - 6.8) <= 1.1 and 4.5 <= em * 10 <= 16.5))
TD = pd.DataFrame(rows)
with pd.option_context("display.float_format", "{:.3g}".format):
print(TD.groupby(["roi_reading", "threshold", "field"]).agg(threshold_vcm=("threshold_vcm", "first"), emission_min=("emission_pct", "min"),
emission_max=("emission_pct", "max"), mV_min=("mV", "min"), mV_max=("mV", "max"), ok=("ok", "sum"), n=("ok", "size")).to_string())
print("with the paper's own threshold, 57.7 V/cm: emission %.2f to %.2f %%" % (TD.sine.min() * 57.7 / 1000, TD.sine.max() * 57.7 / 1000))
print("combinations in range, by slope reading and field reading:", TD[TD.ok].groupby(["slope", "field"]).size().to_dict(),
"of", TD.groupby(["slope", "field"]).size().to_dict())
threshold_vcm emission_min emission_max mV_min mV_max ok n
roi_reading threshold field
aaa t1 amplitude 77.8 0.432 0.629 4.32 6.29 2 8
peak_to_peak 156 0.863 1.26 8.63 12.6 0 8
t5 amplitude 84.1 0.467 0.68 4.67 6.8 4 8
peak_to_peak 168 0.933 1.36 9.33 13.6 0 8
baa t1 amplitude 87.2 0.484 0.705 4.84 7.05 4 8
peak_to_peak 174 0.968 1.41 9.68 14.1 0 8
t5 amplitude 83.7 0.465 0.677 4.65 6.77 4 8
peak_to_peak 167 0.929 1.35 9.29 13.5 0 8
caa t1 amplitude 76 0.422 0.614 4.22 6.14 2 8
peak_to_peak 152 0.844 1.23 8.44 12.3 0 8
t5 amplitude 90.1 0.5 0.728 5 7.28 4 8
peak_to_peak 180 1 1.46 10 14.6 0 8
with the paper's own threshold, 57.7 V/cm: emission 0.32 to 0.47 %
combinations in range, by slope reading and field reading: {('b', 'amplitude'): 20} of {('a', 'amplitude'): 24, ('a', 'peak_to_peak'): 24, ('b', 'amplitude'): 24, ('b', 'peak_to_peak'): 24}
Result: blocked, metadata insufficient. The estimate multiplies the 2-kHz AP threshold in V/cm by the sine sensitivity. The threshold in V/cm cannot be formed from the files: they give labelled generator voltages and "1V translates into 36 V/cm" without saying whether of the amplitude or the peak-to-peak value, and they say the delivered amplitude is "somewhat smaller", the calibration not being deposited. The outcome changes with the field reading. Shown, not deciding: with the signed-regression sine sensitivity (5.55 % per kV/cm) and the threshold under the amplitude reading (76.0 V/cm for 1 AP, 90.1 for 5 APs, crop ROI), the emission change is 0.42 to 0.50 %, 4.2 to 5.0 mV; under the peak-to-peak reading it is 0.84 to 1.46 % over the readings; with the paper's own 57.7 V/cm, 0.32 to 0.47 %. In range are 20 of the 24 combinations with the per-cycle slope under the amplitude reading, and none otherwise. The first round recorded this claim as not reproduced.
The Figure 5 and S3 cells (fig5-sine-asymmetry, fig5-whole-cell-depolarization): not identifiable¶
Both claims are about two named cells (C7 and C8 of Figure 6B). The files do not say which cells these are, and the Figure 5 cell's 12 sine sessions match no cell (five at most, in 2_20240522). Shown here, not deciding anything, are the same measures for the nine cells: the correlation of the two ROIs' oscillations and the LOESS fit at the envelope peak (span 2 ms, scaled). The Pearson correlation of the two ROIs' oscillations takes in everything the LOESS leaves, not only the carrier; the carrier alone is compared in fig6-sine-less-than-pulse (fitted phases).
rows = []
for c in cells6:
w = [x for x in WHOLE if x["cell"] == c and x["span"] == 2.0]
o1, o2 = w[0]["osc"], w[1]["osc"]; m = min(len(o1), len(o2))
rr = stats.pearsonr(o1[:m], o2[:m])
rows.append(dict(cell=c, r_roi1_roi2=rr[0], p=rr[1], loess_peak_roi1=w[0]["loess_at_peak"], loess_peak_roi2=w[1]["loess_at_peak"]))
with pd.option_context("display.float_format", "{:.3g}".format):
print(pd.DataFrame(rows).to_string())
x = SE[(SE.o5 == "a") & (SE.span == 2.0)].pivot(index="cell", columns="roi", values=["pulse", "sine_a"])
x[("pulse_asym_pct", "")] = 100 * (1 - x["pulse"].min(axis=1) / x["pulse"].max(axis=1)); x[("sine_asym_pct", "")] = 100 * (1 - x["sine_a"].min(axis=1) / x["sine_a"].max(axis=1))
print(x.round(1).to_string())
cell r_roi1_roi2 p loess_peak_roi1 loess_peak_roi2
0 2_20240522 0.376 3.34e-18 -0.287 -0.113
1 2_20240523 0.683 6.32e-70 -0.858 -0.435
2 3_20240522 0.663 1.43e-64 -0.174 0.0634
3 3_20240523 0.066 0.141 -0.974 -0.73
4 4_20240522 0.31 1.26e-12 -1.76 -1.5
5 4_20240523 -0.145 0.00118 -0.316 -0.551
6 5_20240523 0.0556 0.215 -0.242 -0.168
7 6_20240523 -0.132 0.00318 0.0508 0.98
8 7_20240523 0.0266 0.554 1.75 1.65
pulse sine_a pulse_asym_pct sine_asym_pct
roi ROI1 ROI2 ROI1 ROI2
cell
2_20240522 14.8 16.6 0.7 8.8 10.8 91.9
2_20240523 27.1 8.3 11.3 0.0 69.2 99.7
3_20240522 21.5 6.9 6.0 3.3 68.1 45.2
3_20240523 16.4 13.9 6.7 4.7 15.4 29.9
4_20240522 20.8 15.9 8.0 6.8 23.4 14.0
4_20240523 7.5 3.1 1.0 1.4 58.1 31.5
5_20240523 19.0 24.2 3.9 4.1 21.6 4.6
6_20240523 29.2 23.4 4.1 12.3 19.9 66.8
7_20240523 26.0 17.7 14.8 2.0 32.0 86.3
fig3b-low-conductivity-halves-threshold, fig3-extrapolated-threshold¶
Blocked: no dandiset holds 250-ms sessions in the low-conductivity (NaGlu) solution (the solution counts in the data check: 001032, 001033 and 001035 are all NaCl; 001030 holds only the 2-s ramped waves). The extrapolation depends on these thresholds.
AP detector variants¶
The claims file's O3 asks for a "width at half height of 2 to 20 ms". The width can be taken at half
the peak's dF (its height above F0) or at half its prominence (its rise above the neighbouring
troughs). Each alone fails on these traces in a way the other catches: at half the dF a small
wiggle on the tail of an AP counts as another AP 5 ms later, and at half the prominence (the first
round's detector) small oscillations on a raised plateau count as APs, so that no session can show
a sustained depolarization. The detector used above asks for 2 to 20 ms under both (tisine.detect_aps).
This cell reruns the time-lapse claims under each single measure, for ROI c and F0 a.
VAR = []
for wr in ["both", "height", "prominence"]:
T = {d: ts.session_table(REC[d], d, wr) for d in ["001032", "001033", "001035", "001271"]}
F = pd.concat([T[d] for d in ["001032", "001033", "001035"]], ignore_index=True)
F = F[F.kind.isin(["sine", "mod"]) & F.carrier.isin([2000.0, 20000.0])].copy()
F["stim"] = np.where(F["kind"] == "mod", "mod", "pure") + "_" + (F.carrier / 1000).astype(int).astype(str) + "k"
for o3 in "ab":
r = ("c", "a", o3); g = rd(F, r)
f1 = fig1(g); th = thresholds(g)
pv = {}
for car, k in [("2k", "t5"), ("20k", "t1"), ("20k", "t5")]:
m = th[th.stim == f"mod_{car}"][k].dropna(); p_ = th[th.stim == f"pure_{car}"][k].dropna()
pv[f"p_mod_{car}_{k}"] = stats.ttest_ind(m, p_).pvalue if len(m) > 1 and len(p_) > 1 else np.nan
c = g.groupby(["neuron", "carrier", "volts", "kind"]).n_ap.mean().unstack("kind").dropna()
bb = block(g, "b"); bt = beats(g, "b", F); tt = ttx(r, "b", T["001271"]); tt = tt[tt.common_volts > 0]
f3 = judge(fig3c("c", "a", o3, "c7_free", wr), "b")
VAR.append(dict(width=wr, o3=o3, excited_2k=f1["excited_2k"], excited_20k=f1["excited_20k"], latency_med=f1["median_latency_at_threshold"],
rho_neg=f"{f1['rho_neg']}/{f1['n_rho']}", max_rate_med=f1["max_rate_median"], block_b=f"{int(bb.sustained.sum())}/{len(bb)}",
mod_max_aps=bt["max_aps"], mod_near_a=bt["near_a"], mod_sustained_b=bt["sustained"], pure_more=f"{int((c.sine > c['mod']).sum())}/{int((c['mod'] > c.sine).sum())}",
ttx_ok=f"{int(tt[tt.group == 'TTX'].ok.sum())}/{int((tt.group == 'TTX').sum())}", vehicle_ok=f"{int(tt[tt.group == 'vehicle'].ok.sum())}/{int((tt.group == 'vehicle').sum())}",
fig3c_block=sum(v["block"] for v in f3.values()), **pv))
with pd.option_context("display.float_format", "{:.3g}".format):
print(pd.DataFrame(VAR).to_string())
print("block_b: neuron-frequencies with sustained depolarization (O4 b) at the highest voltage / tested above their best-rate voltage; "
"pure_more: pairs with more APs for pure / for modulated waves; fig3c_block: neurons of 001030 with the block (O4 b)")
width o3 excited_2k excited_20k latency_med rho_neg max_rate_med block_b mod_max_aps mod_near_a mod_sustained_b pure_more ttx_ok vehicle_ok fig3c_block p_mod_2k_t5 p_mod_20k_t1 p_mod_20k_t5 0 both a 96.4 88.2 60.8 34/37 41.7 0/14 8 0.964 0 44/2 5/6 11/13 0 0.118 0.0337 0.133 1 both b 96.4 94.1 59.8 35/38 40.1 0/14 8 0.964 0 44/2 5/6 11/13 0 0.118 0.102 0.133 2 height a 96.4 88.2 60.8 34/37 136 0/20 16 0.895 1 41/5 5/6 11/13 0 0.369 0.0337 0.272 3 height b 96.4 94.1 59.8 35/38 182 0/22 18 0.873 1 41/5 5/6 11/13 0 0.374 0.102 0.272 4 prominence a 96.4 88.2 66.5 40/40 42.9 0/11 10 0.951 0 47/0 6/6 11/13 0 0.09 0.0205 0.0305 5 prominence b 96.4 94.1 66 41/41 42.3 0/11 10 0.951 0 47/0 6/6 11/13 0 0.09 0.077 0.0225 block_b: neuron-frequencies with sustained depolarization (O4 b) at the highest voltage / tested above their best-rate voltage; pure_more: pairs with more APs for pure / for modulated waves; fig3c_block: neurons of 001030 with the block (O4 b)
Summary¶
| Claim | Outcome |
|---|---|
| ap-selection | partially reproduced (assumed: crop ROI) |
| fig1-pure-sine-excites-all | partially reproduced |
| fig1-excitation-block | not reproduced |
| fig1-20khz-needs-twice-the-field | blocked: metadata insufficient (calibration per frequency) |
| fig2d-5ap-above-1ap | reproduced (assumed: crop ROI) |
| fig2d-modulated-thresholds-higher | partially reproduced (assumed: crop ROI) |
| fig2-modulated-locks-to-beats | reproduced (assumed: envelope phase, O4 level) |
| fig2c-pure-more-aps | reproduced (assumed: crop ROI) |
| fig3b-low-conductivity-halves-threshold | blocked: other (data deposited nowhere) |
| fig3-extrapolated-threshold | blocked: other (depends on fig3b) |
| fig3-not-heating | blocked: methods underspecified |
| fig3c-sustained-firing-2s | partially reproduced |
| figs1-di8-same-as-fluovolt | partially reproduced (assumed: crop ROI) |
| figs2-ttx-blocks-depolarization | partially reproduced (assumed: vehicle split, crop ROI) |
| fig4-pooled-time-constants | blocked: metadata insufficient (no ROIs, experiments not identifiable) |
| fig4-orientation-and-poles | blocked: metadata insufficient |
| fig5-sine-asymmetry | blocked: metadata insufficient |
| fig5-whole-cell-depolarization | blocked: metadata insufficient |
| fig6-sine-less-than-pulse | blocked: metadata insufficient (no ROIs, cells not identifiable, fields not measured) |
| threshold-depolarization-estimate | blocked: metadata insufficient (calibration) |
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 data are the ten dandisets the paper lists, not 001048, which names the paper but holds
another study. Every claim is
data: elsewhere. - The soma ROI: the crop the depositors cut around the soma (16 x 16 pixels) is used as the main reading, with the claims file's Otsu readings shown; the paper's ROIs are not deposited.
- AP detection, sustained depolarization and the windows follow the claims file's readings; the paper identified APs by eye. The AP width at half height is required to be 2 to 20 ms both at half the peak's dF and at half its prominence (an assumption; each alone is shown as a variant).
- Where the ROI readings give different outcomes, the crop ROI decides, as an assumption (protocol step 1.6), where the claims file's header would block.
- Where the O4 readings give different outcomes (fig2-modulated-locks-to-beats only), half the neuron's median AP height (O4 b) decides, as an assumption.
- Fields are labelled generator voltages, not calibrated amplitudes; the conversion to V/cm uses the files' 36 V/cm per volt under two readings.
- Strobe ROIs are 12 x 12 pixel boxes at the extremes of the smoothed ratio image (6 x 6 and 20 x 20 as a variant), flash intensity is scaled to the image median (no debris ROI is marked), and the fields are the paper's 257 and 136 V/cm, scaled by the generator voltage for the two cells run at other voltages.
- The 19 Figure 4 experiments and the Figure 5 and S3 cells cannot be identified in the files.
What was found about the datasets:
- 001048, the dandiset that names the paper, holds HEK-cell TIRF recordings of another study.
- 001034 is a copy of 001032: all 178 recordings are identical, under other subject names.
- 001030 states mouse neurons (C57BL/6), where the paper used rat neurons; 001035 calls its neurons "C57BL/6 rat"; 001030 makes each session a separate subject.
- The voltages are not calibrated: the files say the true amplitude is smaller by frequency-dependent factors "reported separately", which are not in the dandisets.
- No low-conductivity (NaGlu) recordings with 250-ms waves are deposited, so Figures 3A and 3B cannot be tested; 001272 holds 3 of the 8 Di-8-ANEPPS neurons.
- The dandisets hold stimuli the paper does not describe (20-Hz to 10-kHz sine waves in 001032, 001033 and 001035, 40-kHz sine waves in 001033 and 100-kHz ones in 001035; 100-kHz pulses in 001031; TTX strobe sessions in 001274).
- The TTX session descriptions say TTX "caused sustained depolarization", the opposite of the paper and of the dandiset's own description; the vehicle neurons carry no mark of when the vehicle was applied; the Di-8-ANEPPS sessions of 001272 are all described as "vehicle control".
- Strobe files hold no ROIs, no cell shape or orientation labels and no stimulus waveforms.
- 001273 gives the delivered amplitudes of its pulses in a note (25, 62, 79 and 102 V for 40, 80, 100 and 130 V), where the paper gives 140 V/cm for every experiment.