Reproducing Quattrocolo & Maccaferri (2013) from CRCNS hc-9¶
Paper. Quattrocolo G, Maccaferri G (2013). Novel GABAergic circuits mediating
excitation/inhibition of Cajal-Retzius cells in the developing hippocampus. J Neurosci
33(13):5486-5498. doi:10.1523/jneurosci.5680-12.2013.
The claims file (claims/crcns-hc-9/) was written from the version of record, which has no
supplement and links no code.
Data set. CRCNS hc-9, "Firing pattern of O-LM cells in mouse hippocampal CA1"
(doi:10.6080/K0HQ3WTH), read from the public bucket at
https://crcnsarchive.s3-us-west-2.amazonaws.com/hc-9/live/. CRCNS data sets are not versioned,
so the version is recorded as unversioned. The data are pinned by MD5: the one data file,
Maccaferri_data.zip, has MD5 4e3175055681b99ad5be06ba542a93b8, as the data set's
checksums.md5 gives it, and the data check below confirms the MD5 of the zip as read. This is the
only published copy of the data (the bucket's hc-9/old/ folder holds only earlier documentation),
so there is no version to choose between; the paper cites no data set. The description read is
version 0.7 (Dec 2, 2016) of live/docs/crcns_hc-9_data_description.pdf.
What the data set holds. Current-clamp responses of three O-LM interneurons to 1 s current steps (15 sweeps each), the firing patterns of the paper's Figure 7. None of the recordings from Cajal-Retzius cells behind the paper's main findings (Figures 1 to 6 and 8 to 12) are in it, and the light-evoked responses of Figure 7C and 7D are not either. The zip is opened in place over range requests; nothing is written to disk.
Claims. One of the 26 claims can be tested.
| Claim | Tested or blocked |
|---|---|
| fig7a-olm-sag | tested |
| fig1-spontaneous-amplitude | blocked, data missing |
| fig4-upsc-amplitude | blocked, data missing |
| fig4-upsc-kinetics | blocked, data missing |
| paired-presynaptic-neurogliaform | blocked, data missing |
| upsc-gabazine-block | blocked, data missing |
| fig5-no-gabab-component | blocked, data missing |
| fig6-spontaneous-larger-than-unitary | blocked, data missing |
| fig6-spontaneous-faster-than-unitary | blocked, data missing |
| fig7-chr2-olm-light-firing | blocked, data missing |
| fig8-light-evoked-psc | blocked, data missing |
| fig9-opsc-rise-time | blocked, data missing |
| fig9-opsc-amplitude-decay | blocked, data missing |
| fig10-dhpg-increases-frequency | blocked, data missing |
| fig10-dhpg-amplitude-unchanged | blocked, data missing |
| fig10-reduced-slice-no-effect | blocked, data missing |
| fig11-slm-inhibition | blocked, data missing |
| fig11-slm-inhibition-gabazine | blocked, data missing |
| fig11-slm-excitation | blocked, data missing |
| fig11d-baseline-frequency-dichotomy | blocked, data missing |
| gabazine-baseline-firing | blocked, data missing |
| fig12a-slow-conductance-depolarized | blocked, data missing |
| fig12a-slow-conductance-hyperpolarized | blocked, data missing |
| fig12b-fast-conductance-depolarized | blocked, data missing |
| fig12b-fast-conductance-hyperpolarized | blocked, data missing |
| no-age-dependence | blocked, data missing |
import io
import itertools
import time
import zipfile
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
from scipy import stats
import hc9
pd.set_option("display.width", 200)
pd.set_option("display.max_columns", 40)
started = time.time()
cells = hc9.load_cells()
print(f"loaded {len(cells)} cells from CRCNS {hc9.DATASET} ({hc9.VERSION}), {hc9.ZIP_URL}, "
f"in {time.time() - started:.1f} s")
loaded 3 cells from CRCNS hc-9 (unversioned), https://crcnsarchive.s3-us-west-2.amazonaws.com/hc-9/live/Maccaferri_data.zip, in 3.3 s
Data check¶
The data set is one zip in the CRCNS bucket. Its MD5 is compared with the data set's own
checksums.md5, and the zip's members are listed. Each check of the description against the
values is stated under its output.
published = hc9.published_checksums()
print("checksums.md5:")
for k, v in published.items():
print(f" {v} {k}")
md5 = hc9.zip_md5()
print(f"\nMD5 of Maccaferri_data.zip as read now: {md5}")
print("matches checksums.md5:", md5 == published["Maccaferri_data.zip"],
"; matches the MD5 pinned in hc9.py:", md5 == hc9.EXPECTED_ZIP_MD5)
print("\nMembers of the zip:")
print(hc9.zip_listing().to_string(index=False))
print("\nMD5 of each ATF member:")
for name, c in cells.items():
print(f" {c['raw_md5']} {name}_raw.atf")
checksums.md5: 32366f0262dc21cb97cb1f1efaade1a7 docs/crcns_hc-9_data_description.docx 7b36e1831089cad1aa854b0e1bf9b011 docs/crcns_hc-9_data_description.pdf 23542905cd6ea65f31f492f63fae0bfc filelist.txt 4e3175055681b99ad5be06ba542a93b8 Maccaferri_data.zip
MD5 of Maccaferri_data.zip as read now: 4e3175055681b99ad5be06ba542a93b8 matches checksums.md5: True ; matches the MD5 pinned in hc9.py: True Members of the zip:
member bytes date
Maccaferri_data/ 0 2016-11-23
Maccaferri_data/1May2012_P3_raw.atf 7388537 2014-05-06
Maccaferri_data/20Sept2011_P2_raw.atf 7388498 2014-05-06
Maccaferri_data/24October2012_C2_raw.atf 7385994 2014-05-06
MD5 of each ATF member:
3a1a661c64058505f0d221c3ba309b6b 1May2012_P3_raw.atf
c4097c125efe088d6666fad6371d7947 20Sept2011_P2_raw.atf
a4601337699cf96ea8aa3b0413a15aa4 24October2012_C2_raw.atf
for name, c in cells.items():
atf = c["atf"]
print(f"== {name}: {atf['version']}, {atf['n_header']} optional header records, {atf['n_cols']} columns")
for h in atf["header"]:
print(" ", h[:160])
print(" columns:", atf["columns"][:5], "...", atf["columns"][-2:])
dt = np.diff(c["t"])
print(f" samples {len(c['t'])}, time {c['t'][0]:.4f} to {c['t'][-1]:.4f} s, "
f"interval {dt.min()*1e6:.3f} to {dt.max()*1e6:.3f} us, sampling {c['fs']:.0f} Hz")
print(f" sweeps {c['v'].shape[1]}; voltage {c['v'].min():.1f} to {c['v'].max():.1f} mV; "
f"current {c['i'].min():.1f} to {c['i'].max():.1f} pA")
== 1May2012_P3: ATF 1.0, 8 optional header records, 31 columns
"AcquisitionMode=Episodic Stimulation"
"Comment="
"YTop=200,4291.76"
"YBottom=-200,-4334.72"
"SyncTimeUnits=12.5"
"SweepStartTimesMS=0.000,5000.000,10000.000,15000.000,20000.000,25000.000,30000.000,35000.000,40000.000,45000.000,50000.000,55000.000,60000.000,65000.000,70000.
"SignalsExported=2CC mV,AO #1"
"Signals=" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO
columns: ['Time (s)', 'Trace #1 (mV)', 'Trace #1 (pA)', 'Trace #2 (mV)', 'Trace #2 (pA)'] ... ['Trace #15 (mV)', 'Trace #15 (pA)']
samples 39999, time 0.0000 to 1.9999 s, interval 50.000 to 50.000 us, sampling 20000 Hz
sweeps 15; voltage -102.8 to 33.9 mV; current -100.0 to 250.0 pA
== 20Sept2011_P2: ATF 1.0, 8 optional header records, 31 columns
"AcquisitionMode=Episodic Stimulation"
"Comment="
"YTop=200,4291.76"
"YBottom=-200,-4334.72"
"SyncTimeUnits=12.5"
"SweepStartTimesMS=0.000,5000.000,10000.000,15000.000,20000.000,25000.000,30000.000,35000.000,40000.000,45000.000,50000.000,55000.000,60000.000,65000.000,70000.
"SignalsExported=2CC mV,AO #1"
"Signals=" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO
columns: ['Time (s)', 'Trace #1 (mV)', 'Trace #1 (pA)', 'Trace #2 (mV)', 'Trace #2 (pA)'] ... ['Trace #15 (mV)', 'Trace #15 (pA)']
samples 39999, time 0.0000 to 1.9999 s, interval 50.000 to 50.000 us, sampling 20000 Hz
sweeps 15; voltage -99.0 to 29.9 mV; current -100.0 to 250.0 pA
== 24October2012_C2: ATF 1.0, 8 optional header records, 31 columns
"AcquisitionMode=Episodic Stimulation"
"Comment="
"YTop=200,4295.29"
"YBottom=-200,-4340.95"
"SyncTimeUnits=12.5"
"SweepStartTimesMS=0.000,5000.000,10000.000,15000.000,20000.000,25000.000,30000.000,35000.000,40000.000,45000.000,50000.000,55000.000,60000.000,65000.000,70000.
"SignalsExported=2CC mV,AO #1"
"Signals=" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO #1" "2CC mV" "AO
columns: ['Time (s)', 'Trace #1 (mV)', 'Trace #1 (pA)', 'Trace #2 (mV)', 'Trace #2 (pA)'] ... ['Trace #15 (mV)', 'Trace #15 (pA)']
samples 39999, time 0.0000 to 1.9999 s, interval 50.000 to 50.000 us, sampling 20000 Hz
sweeps 15; voltage -104.6 to 38.3 mV; current -100.0 to 250.0 pA
The ATF header agrees with the description: line 2 gives the number of further header records (8) and of columns (31), the 11th line names the columns, time is in seconds, voltage in mV and the injected current in pA. Each file holds 15 sweeps of 2 s sampled at 20 kHz (the Methods give 10 to 20 kHz), with sweep starts 5 s apart (the description says 3 to 5 s).
steps = pd.concat({name: hc9.step_table(c) for name, c in cells.items()}, names=["cell"]).reset_index(level=0)
print("Step level, onset and offset of every sweep, read from the current channel:")
print(steps.pivot(index="sweep", columns="cell", values="step_pA").to_string())
print("\nOnsets (s):", sorted(steps.on_s.dropna().round(5).unique()),
" offsets (s):", sorted(steps.off_s.dropna().round(5).unique()),
" holding current (pA):", sorted(steps.holding_pA.round(2).unique()))
levels = sorted(steps.step_pA.unique())
print("\nLevels present:", levels)
print("Increment:", np.unique(np.diff(levels)), "pA")
stated = list(range(-100, 251, 30))
print("Description, -100:30:+250 pA, would give", len(stated), "levels:", stated)
Step level, onset and offset of every sweep, read from the current channel: cell 1May2012_P3 20Sept2011_P2 24October2012_C2 sweep 1 -100.0 -100.0 -100.0 2 -75.0 -75.0 -75.0 3 -50.0 -50.0 -50.0 4 -25.0 -25.0 -25.0 5 0.0 0.0 0.0 6 25.0 25.0 25.0 7 50.0 50.0 50.0 8 75.0 75.0 75.0 9 100.0 100.0 100.0 10 125.0 125.0 125.0 11 150.0 150.0 150.0 12 175.0 175.0 175.0 13 200.0 200.0 200.0 14 225.0 225.0 225.0 15 250.0 250.0 250.0 Onsets (s): [np.float64(0.07505)] offsets (s): [np.float64(1.07505)] holding current (pA): [np.int64(0)] Levels present: [np.float64(-100.0), np.float64(-75.0), np.float64(-50.0), np.float64(-25.0), np.float64(0.0), np.float64(25.0), np.float64(50.0), np.float64(75.0), np.float64(100.0), np.float64(125.0), np.float64(150.0), np.float64(175.0), np.float64(200.0), np.float64(225.0), np.float64(250.0)] Increment: [25.] pA Description, -100:30:+250 pA, would give 12 levels: [-100, -70, -40, -10, 20, 50, 80, 110, 140, 170, 200, 230]
The description gives the current steps as -100:30:+250 pA. The current channel shows 15 steps from -100 to +250 pA in 25 pA increments, one of them 0 pA, in all three files; -100:30:+250 would give 12 levels ending at +230 pA. The data are trusted here: the description is wrong on the increment, and this is recorded as a finding. Each step lasts 1 s, from 75 ms to 1.075 s; the description says the injection "does not usually start until about half a second into the recording", which does not hold for these files. Both the -100 pA and the +50 pA steps of the Figure 7 legends are present in every file.
The description names the three cells but not their animals. The superseded description kept
in the same bucket (hc-9/old/) gives them:
old = hc9.old_description_text()
lines = [l for l in old.splitlines() if "_raw" in l]
for l in lines:
print(l)
print("\nFigure 7 legend (paper): A, two O-LM cells, P13 and P14, from CXCR4-EGFP mice; B to D, an O-LM "
"cell, P15, from the Sst-IRES-Cre;ChR2(H134R)-EYFP mouse.")
1May2012_P3_raw comes from a P13 CXCR4-EGFP mouse 20Sept2011_P2_raw comes P14 CXCR4-EGFP mouse 24October2012_C2_raw comes from a P15 Sst-IRES-Cre; ChR2(H134R)-EYFP mouse - May2012_P3_raw.atf comes from a P13 CXCR4-EGFP mouse - 20Sept2011_P2_raw.atf comes P14 CXCR4-EGFP mouse - 24October2012_C2_raw.atf comes from a P15 Sst-IRES-Cre; ChR2(H134R)-EYFP mouse Figure 7 legend (paper): A, two O-LM cells, P13 and P14, from CXCR4-EGFP mice; B to D, an O-LM cell, P15, from the Sst-IRES-Cre;ChR2(H134R)-EYFP mouse.
The ages and mouse lines of the old description match the Figure 7 legends: 1May2012_P3 (P13) and 20Sept2011_P2 (P14) are consistent with the two cells of Figure 7A, and 24October2012_C2 (P15, Sst-IRES-Cre;ChR2) with the cell of Figure 7B to D. This is consistent with the data set holding the three Figure 7 cells, as its description says, but the files themselves carry no age, line or panel. The old description also calls the files ABF; they are ATF text files.
Next, a check for recordings deposited more than once: every pair of sweeps, within and across files, is compared.
rows = []
keys = [(n, k) for n, c in cells.items() for k in range(c["v"].shape[1])]
for (n1, k1), (n2, k2) in itertools.combinations(keys, 2):
a, b = cells[n1]["v"][:, k1], cells[n2]["v"][:, k2]
rows.append({"a": f"{n1}#{k1+1}", "b": f"{n2}#{k2+1}", "same_cell": n1 == n2,
"identical": bool(np.array_equal(a, b)), "max_abs_diff_mV": float(np.abs(a - b).max()),
"r": float(np.corrcoef(a, b)[0, 1])})
pairs = pd.DataFrame(rows)
print(f"{len(pairs)} pairs of sweeps; identical pairs: {int(pairs.identical.sum())}")
print("Smallest maximum difference between two sweeps (mV):", round(pairs.max_abs_diff_mV.min(), 3))
print("\nAcross cells, the same step level: correlation of the voltage traces")
same_level = []
for lvl in levels:
ks = {n: hc9.sweep_for(c, lvl) for n, c in cells.items()}
for n1, n2 in itertools.combinations(cells, 2):
a, b = cells[n1]["v"][:, ks[n1]], cells[n2]["v"][:, ks[n2]]
same_level.append({"level_pA": lvl, "pair": f"{n1} / {n2}", "r": np.corrcoef(a, b)[0, 1],
"max_abs_diff_mV": np.abs(a - b).max()})
same_level = pd.DataFrame(same_level)
print(same_level.groupby("pair")[["r", "max_abs_diff_mV"]].agg(["min", "max"]).round(3).to_string())
990 pairs of sweeps; identical pairs: 0
Smallest maximum difference between two sweeps (mV): 2.704
Across cells, the same step level: correlation of the voltage traces
r max_abs_diff_mV
min max min max
pair
1May2012_P3 / 20Sept2011_P2 -0.006 0.990 3.058 96.112
1May2012_P3 / 24October2012_C2 0.126 0.993 3.296 104.083
20Sept2011_P2 / 24October2012_C2 0.230 0.993 2.704 101.483
No two sweeps are identical, within a file or across files, and the three files differ at every step level, so no recording is deposited twice. Each cell is one recording from one animal, so results below are given per cell.
The voltage must share the clock of the current channel. Two checks: the delay from the step onset in the current channel to the first departure of the voltage from its baseline, and the voltage change in the step window against the step level, with the window as the current channel gives it and shifted in time by 1.1 s (after the step).
rows = []
for name, c in cells.items():
t, v = c["t"], c["v"]
on, off = hc9.step_window(c)
st = hc9.step_table(c)
pre = t < on
for k, lvl in enumerate(st.step_pA):
if lvl == 0:
continue
base, sd = v[pre, k].mean(), v[pre, k].std(ddof=1)
dev = np.flatnonzero((np.abs(v[:, k] - base) > 5 * sd) & (t >= on - 0.02))
rows.append({"cell": name, "level_pA": lvl, "voltage_departs_ms_after_onset": (t[dev[0]] - on) * 1e3})
lag = pd.DataFrame(rows)
print("Delay from the current step's onset to the voltage leaving baseline by 5 SD (ms):")
print(lag.groupby("cell").voltage_departs_ms_after_onset.describe()[["count", "min", "50%", "max"]].round(2).to_string())
def window_change(c, shift):
on, off = hc9.step_window(c)
t, v = c["t"], c["v"]
a, b = on + shift, min(off + shift, t[-1])
w = (t >= a + 0.2 * (b - a)) & (t < b)
pre = t < on
# median over the window, so spikes on depolarizing steps do not dominate
return np.median(v[w], axis=0) - v[pre].mean(axis=0)
rows = []
for name, c in cells.items():
lv = hc9.step_table(c).step_pA.to_numpy()
neg = lv < 0
for label, shift in [("step window", 0.0), ("shifted by 1.1 s", 1.1 - 0.075)]:
dv = window_change(c, shift)
rho, p = stats.spearmanr(lv[neg], dv[neg])
rows.append({"cell": name, "window": label, "spearman_rho_neg_steps": rho, "p": p,
"dV_at_-100pA_mV": dv[lv == -100][0], "dV_at_+250pA_mV": dv[lv == 250][0]})
print("\nVoltage change (median over the window's last 80%, minus baseline) against step level:")
print(pd.DataFrame(rows).round(4).to_string(index=False))
Delay from the current step's onset to the voltage leaving baseline by 5 SD (ms):
count min 50% max
cell
1May2012_P3 14.0 0.25 1.52 4.00
20Sept2011_P2 14.0 0.35 0.55 3.35
24October2012_C2 14.0 0.10 0.75 2.55
Voltage change (median over the window's last 80%, minus baseline) against step level:
cell window spearman_rho_neg_steps p dV_at_-100pA_mV dV_at_+250pA_mV
1May2012_P3 step window 1.0 0.0 -24.4557 17.9590
1May2012_P3 shifted by 1.1 s -0.2 0.8 1.9359 -3.1591
20Sept2011_P2 step window 1.0 0.0 -28.0433 17.4155
20Sept2011_P2 shifted by 1.1 s -0.2 0.8 0.3807 -3.3608
24October2012_C2 step window 1.0 0.0 -34.2732 24.1759
24October2012_C2 shifted by 1.1 s -0.8 0.2 4.9114 -5.0050
The voltage leaves its baseline 0.10 to 4.00 ms after the current step's onset in every sweep (medians 0.55 to 1.52 ms per cell), and in the step window the voltage change follows the hyperpolarizing step level monotonically (Spearman rho 1.0 in each cell), while in the window shifted past the step it does not (rho -0.8 to -0.2, p 0.2 to 0.8). The voltage and current channels share one clock. The three cells' pre-step potentials, about -63 mV, are given in the next section.
fig7a-olm-sag¶
O-LM interneurons recorded in CXCR4-EGFP slices show a sag in their voltage response to hyperpolarizing current pulses (pulses of -100 and +50 pA in Figure 7A), and the ChR2-expressing O-LM cell of Figure 7D has a similar firing pattern.
Criterion (claims file): the -100 pA step, judged under two readings. (a) The relaxation from the
peak hyperpolarization during the step to the mean voltage over its last 100 ms is larger than
twice the SD of the pre-step baseline. (b) The sag ratio (V_end - V_min) / (V_baseline - V_min) is
at least 0.1. Every cell, the Figure 7D cell included, is judged alike. Reproduced when every cell
shows a sag, partially when most do. The parameters are in hc9.py. Two details are standard choices: the pre-step baseline
is all of the sweep before the step (0 to 75 ms), and V_min is the minimum of a 1 ms running
mean during the step (the raw minimum is a variant; the noise checks below show why it is not
the main measure).
on_off = {name: hc9.step_window(c) for name, c in cells.items()}
panel = {"1May2012_P3": "7A (P13, CXCR4-EGFP)", "20Sept2011_P2": "7A (P14, CXCR4-EGFP)",
"24October2012_C2": "7D (P15, Sst-Cre;ChR2)"}
rows = []
for name, c in cells.items():
on, off = on_off[name]
k = hc9.sweep_for(c, hc9.HYPER_STEP_PA)
m = hc9.sag_measures(c["t"], c["v"][:, k], on, off)
rows.append({"cell": name, "panel (old description)": panel[name], **m})
sag = pd.DataFrame(rows)
sag["a_threshold_mV"] = hc9.SD_FACTOR * sag.sd_base
sag["reading_a"] = sag.sag_mV > sag.a_threshold_mV
sag["reading_b"] = sag.sag_ratio >= hc9.RATIO_THRESHOLD
print("-100 pA step, main measures (V in mV, t_min in ms from step onset):")
print(sag.round(3).to_string(index=False))
print(f"\nReading (a): {int(sag.reading_a.sum())} of {len(sag)} cells show a sag; "
f"reading (b): {int(sag.reading_b.sum())} of {len(sag)}.")
print("Sag over the reading (a) threshold, per cell:", (sag.sag_mV / sag.a_threshold_mV).round(1).tolist())
print("Sag ratio over the reading (b) threshold, per cell:", (sag.sag_ratio / hc9.RATIO_THRESHOLD).round(2).tolist())
-100 pA step, main measures (V in mV, t_min in ms from step onset):
cell panel (old description) v_base sd_base v_min t_min_ms v_end sag_mV sag_over_sd sag_ratio a_threshold_mV reading_a reading_b
1May2012_P3 7A (P13, CXCR4-EGFP) -63.575 0.279 -102.735 116.40 -85.196 17.539 62.909 0.448 0.558 True True
20Sept2011_P2 7A (P14, CXCR4-EGFP) -63.509 0.182 -98.974 88.95 -91.038 7.936 43.674 0.224 0.363 True True
24October2012_C2 7D (P15, Sst-Cre;ChR2) -62.974 0.093 -104.546 125.65 -95.880 8.666 93.160 0.208 0.186 True True
Reading (a): 3 of 3 cells show a sag; reading (b): 3 of 3.
Sag over the reading (a) threshold, per cell: [31.5, 21.8, 46.6]
Sag ratio over the reading (b) threshold, per cell: [4.48, 2.24, 2.08]
How the result changes with the choices the criterion leaves to a standard: V_min from the raw trace or a 1 ms running mean, the end window (50, 100 and 200 ms), the multiple of the baseline SD (2, 3, 5 and 10) and the ratio threshold (0.05 to 0.3).
rows = []
for name, c in cells.items():
on, off = on_off[name]
k = hc9.sweep_for(c, hc9.HYPER_STEP_PA)
for smooth in (0.0, hc9.SMOOTH_S):
for end in (0.05, 0.1, 0.2):
m = hc9.sag_measures(c["t"], c["v"][:, k], on, off, end_window=end, smooth=smooth)
rows.append({"cell": name, "smooth_ms": smooth * 1e3, "end_window_ms": end * 1e3,
"sag_mV": m["sag_mV"], "sag_over_sd": m["sag_over_sd"], "sag_ratio": m["sag_ratio"]})
var = pd.DataFrame(rows)
print(var.round(3).to_string(index=False))
print("\nSmallest sag over baseline SD in any variant:", round(var.sag_over_sd.min(), 1),
"; smallest sag ratio:", round(var.sag_ratio.min(), 3))
for f in (2, 3, 5, 10):
print(f" reading (a) with {f} SD: cells with a sag in every variant:",
int(var.groupby("cell").sag_over_sd.min().gt(f).sum()), "of 3")
for r in (0.05, 0.1, 0.15, 0.2, 0.25, 0.3):
print(f" reading (b) with ratio >= {r}: cells with a sag in every variant:",
int(var.groupby("cell").sag_ratio.min().ge(r).sum()), "of 3")
cell smooth_ms end_window_ms sag_mV sag_over_sd sag_ratio
1May2012_P3 0.0 50.0 17.874 64.110 0.455
1May2012_P3 0.0 100.0 17.636 63.257 0.449
1May2012_P3 0.0 200.0 17.126 61.428 0.436
1May2012_P3 1.0 50.0 17.777 63.762 0.454
1May2012_P3 1.0 100.0 17.539 62.909 0.448
1May2012_P3 1.0 200.0 17.029 61.080 0.435
20Sept2011_P2 0.0 50.0 7.979 43.907 0.225
20Sept2011_P2 0.0 100.0 7.992 43.979 0.225
20Sept2011_P2 0.0 200.0 7.759 42.697 0.218
20Sept2011_P2 1.0 50.0 7.923 43.602 0.223
20Sept2011_P2 1.0 100.0 7.936 43.674 0.224
20Sept2011_P2 1.0 200.0 7.703 42.392 0.217
24October2012_C2 0.0 50.0 8.707 93.599 0.209
24October2012_C2 0.0 100.0 8.765 94.223 0.210
24October2012_C2 0.0 200.0 8.394 90.236 0.201
24October2012_C2 1.0 50.0 8.608 92.535 0.207
24October2012_C2 1.0 100.0 8.666 93.160 0.208
24October2012_C2 1.0 200.0 8.295 89.173 0.200
Smallest sag over baseline SD in any variant: 42.4 ; smallest sag ratio: 0.2
reading (a) with 2 SD: cells with a sag in every variant: 3 of 3
reading (a) with 3 SD: cells with a sag in every variant: 3 of 3
reading (a) with 5 SD: cells with a sag in every variant: 3 of 3
reading (a) with 10 SD: cells with a sag in every variant: 3 of 3
reading (b) with ratio >= 0.05: cells with a sag in every variant: 3 of 3
reading (b) with ratio >= 0.1: cells with a sag in every variant: 3 of 3
reading (b) with ratio >= 0.15: cells with a sag in every variant: 3 of 3
reading (b) with ratio >= 0.2: cells with a sag in every variant: 2 of 3
reading (b) with ratio >= 0.25: cells with a sag in every variant: 1 of 3
reading (b) with ratio >= 0.3: cells with a sag in every variant: 1 of 3
The same measures on the other hyperpolarizing steps, shown as a variant (the legends name the -100 pA step), and the firing on the +50 pA step, the legends' depolarizing pulse. The firing is described only: the criterion judges the similarity of the Figure 7D cell by the sag.
rows = []
for name, c in cells.items():
on, off = on_off[name]
for lvl in (-100, -75, -50, -25):
k = hc9.sweep_for(c, lvl)
m = hc9.sag_measures(c["t"], c["v"][:, k], on, off)
rows.append({"cell": name, "level_pA": lvl, "v_min": m["v_min"], "sag_mV": m["sag_mV"],
"sag_over_sd": m["sag_over_sd"], "sag_ratio": m["sag_ratio"]})
print(pd.DataFrame(rows).round(3).to_string(index=False))
rows = []
for name, c in cells.items():
on, off = on_off[name]
for lvl in (hc9.DEPOL_STEP_PA, 100.0, 250.0):
k = hc9.sweep_for(c, lvl)
sp = hc9.spike_times(c["t"], c["v"][:, k])
inside = sp[(sp >= on) & (sp < off)]
isi = np.diff(inside) * 1e3
rows.append({"cell": name, "level_pA": lvl, "spikes_in_step": len(inside),
"rate_Hz": len(inside) / (off - on),
"first_spike_ms": (inside[0] - on) * 1e3 if len(inside) else np.nan,
"first_isi_ms": isi[0] if len(isi) else np.nan,
"last_isi_ms": isi[-1] if len(isi) else np.nan,
"spikes_outside_step": len(sp) - len(inside)})
fire = pd.DataFrame(rows)
print("\nFiring on depolarizing steps (spikes are upward crossings of 0 mV):")
print(fire.round(1).to_string(index=False))
cell level_pA v_min sag_mV sag_over_sd sag_ratio
1May2012_P3 -100 -102.735 17.539 62.909 0.448
1May2012_P3 -75 -91.979 10.283 76.014 0.360
1May2012_P3 -50 -84.011 5.646 37.544 0.267
1May2012_P3 -25 -74.898 2.747 11.928 0.210
20Sept2011_P2 -100 -98.974 7.936 43.674 0.224
20Sept2011_P2 -75 -89.160 4.256 7.599 0.171
20Sept2011_P2 -50 -81.042 3.176 4.006 0.174
20Sept2011_P2 -25 -72.935 0.674 1.518 0.070
24October2012_C2 -100 -104.546 8.666 93.160 0.208
24October2012_C2 -75 -93.756 4.444 12.007 0.139
24October2012_C2 -50 -85.131 5.334 65.609 0.230
24October2012_C2 -25 -74.939 1.475 9.462 0.117
Firing on depolarizing steps (spikes are upward crossings of 0 mV):
cell level_pA spikes_in_step rate_Hz first_spike_ms first_isi_ms last_isi_ms spikes_outside_step
1May2012_P3 50.0 15 15.0 33.4 56.6 74.7 0
1May2012_P3 100.0 27 27.0 18.2 31.0 43.8 0
1May2012_P3 250.0 48 48.0 7.9 15.6 22.9 0
20Sept2011_P2 50.0 15 15.0 19.2 58.5 93.8 0
20Sept2011_P2 100.0 35 35.0 6.8 23.4 31.6 0
20Sept2011_P2 250.0 83 83.0 3.0 11.0 10.2 0
24October2012_C2 50.0 18 18.0 27.6 49.3 61.5 0
24October2012_C2 100.0 36 36.0 13.7 24.4 29.6 0
24October2012_C2 250.0 60 60.0 5.2 10.6 21.9 0
How large a "sag" the measures find where there is none. Two checks of reading (a)'s false-positive level, with V_min from the raw trace and from the 1 ms running mean. First, each cell's 0 pA sweep, measured in the same step window: no current is injected, so any "sag" is noise. Second, a synthetic trace with no sag: a 40 mV exponential hyperpolarization (time constant 20 ms) in the same window, plus Gaussian noise of 0.2 mV SD, 100 draws from a fixed seed. The 0 pA sweep also gives a reference for the -100 pA sag in mV: the same measure with no current injected, and the SD of the voltage over the 1 s step window, which includes slow fluctuations that the 75 ms pre-step baseline does not. For reading (b) the 0 pA "sag" is also divided by the -100 pA hyperpolarization of the same cell (V_base minus V_min): this is the part of the -100 pA sag ratio that noise of the 0 pA size could make. The sag ratio of the 0 pA sweep itself is printed too, but its denominator, the 0 pA sweep's deflection below baseline, is noise as well, so it does not measure a false-positive level of reading (b).
rows = []
for name, c in cells.items():
on, off = on_off[name]
k = hc9.sweep_for(c, 0.0)
for smooth in (0.0, hc9.SMOOTH_S):
m = hc9.sag_measures(c["t"], c["v"][:, k], on, off, smooth=smooth)
rows.append({"cell": name, "sweep": k + 1, "smooth_ms": smooth * 1e3, "sd_base": m["sd_base"],
"sag_mV": m["sag_mV"], "sag_over_sd": m["sag_over_sd"], "sag_ratio": m["sag_ratio"]})
zero = pd.DataFrame(rows)
print("0 pA sweep, step window as for the -100 pA step:")
print(zero.round(3).to_string(index=False))
rng = np.random.default_rng(20130327)
c0 = cells[hc9.CELLS[0]]
t = c0["t"]
on, off = on_off[hc9.CELLS[0]]
clean = np.where(t >= on, -40.0 * (1 - np.exp(-np.clip(t - on, 0, None) / 0.020)), 0.0) - 63.0
clean[t >= off] = -63.0
rows = []
for draw in range(100):
v = clean + rng.normal(0, 0.2, len(t))
for smooth in (0.0, hc9.SMOOTH_S):
m = hc9.sag_measures(t, v, on, off, smooth=smooth)
rows.append({"draw": draw, "smooth_ms": smooth * 1e3, "sag_over_sd": m["sag_over_sd"],
"sag_ratio": m["sag_ratio"]})
syn = pd.DataFrame(rows)
print("\nSynthetic trace without sag, 100 draws: sag over baseline SD and sag ratio")
print(syn.groupby("smooth_ms")[["sag_over_sd", "sag_ratio"]].describe(percentiles=[0.5])
.round(3).to_string())
print("Fraction of draws passing reading (a) at 2 SD:",
syn.groupby("smooth_ms").sag_over_sd.apply(lambda x: float((x > hc9.SD_FACTOR).mean())).to_dict(),
"; passing reading (b):",
syn.groupby("smooth_ms").sag_ratio.apply(lambda x: float((x >= hc9.RATIO_THRESHOLD).mean())).to_dict())
rows = []
for name, c in cells.items():
on, off = on_off[name]
k = hc9.sweep_for(c, 0.0)
t, v = c["t"], c["v"][:, k]
during = (t >= on) & (t < off)
z = zero[(zero.cell == name) & (zero.smooth_ms == 1.0)].iloc[0]
real = sag.set_index("cell").loc[name]
rows.append({"cell": name, "spikes_0pA": len(hc9.spike_times(t, v)),
"sd_0pA_in_step_mV": v[during].std(ddof=1), "sag_0pA_mV": z.sag_mV,
"sag_0pA_over_sd": z.sag_over_sd, "sag_-100pA_mV": real.sag_mV,
"sag_-100pA_over_sag_0pA": real.sag_mV / z.sag_mV,
"sag_-100pA_over_2sd_0pA_in_step": real.sag_mV / (2 * v[during].std(ddof=1)),
"sag_ratio_0pA_own": z.sag_ratio,
"sag_0pA_over_hyper_-100pA": z.sag_mV / (real.v_base - real.v_min),
"sag_ratio_-100pA": real.sag_ratio})
ref = pd.DataFrame(rows)
print("\nThe -100 pA sag against the same measure on the 0 pA sweep (V_min from the 1 ms running mean):")
print(ref.round(3).to_string(index=False))
print("\nReading (b) on the 0 pA sweeps. sag_ratio_0pA_own divides the 0 pA 'sag' by the 0 pA sweep's own")
print("deflection below baseline, which is noise too, so it is not a sag ratio of a hyperpolarization;")
print("sag_0pA_over_hyper_-100pA divides it by the -100 pA hyperpolarization, the ratio noise alone would")
print("add to the -100 pA sag ratio. Threshold of reading (b):", hc9.RATIO_THRESHOLD)
print(" 0 pA sweeps passing reading (b) on their own deflection:", int((ref.sag_ratio_0pA_own >= hc9.RATIO_THRESHOLD).sum()), "of 3")
print(" noise ratio against the -100 pA hyperpolarization:", ref["sag_0pA_over_hyper_-100pA"].round(3).tolist(),
"; passing:", int((ref["sag_0pA_over_hyper_-100pA"] >= hc9.RATIO_THRESHOLD).sum()), "of 3")
0 pA sweep, step window as for the -100 pA step:
cell sweep smooth_ms sd_base sag_mV sag_over_sd sag_ratio
1May2012_P3 5 0.0 0.286 2.498 8.732 2.885
1May2012_P3 5 1.0 0.286 2.393 8.364 3.146
20Sept2011_P2 5 0.0 0.473 1.425 3.012 0.719
20Sept2011_P2 5 1.0 0.473 1.348 2.848 0.707
24October2012_C2 5 0.0 0.095 2.023 21.258 4.092
24October2012_C2 5 1.0 0.095 1.916 20.136 4.943
Synthetic trace without sag, 100 draws: sag over baseline SD and sag ratio
sag_over_sd sag_ratio
count mean std min 50% max count mean std min 50% max
smooth_ms
0.0 100.0 3.984 0.312 3.438 3.960 4.912 100.0 0.020 0.001 0.017 0.019 0.024
1.0 100.0 0.856 0.091 0.692 0.844 1.274 100.0 0.004 0.000 0.003 0.004 0.006
Fraction of draws passing reading (a) at 2 SD: {0.0: 1.0, 1.0: 0.0} ; passing reading (b): {0.0: 0.0, 1.0: 0.0}
The -100 pA sag against the same measure on the 0 pA sweep (V_min from the 1 ms running mean):
cell spikes_0pA sd_0pA_in_step_mV sag_0pA_mV sag_0pA_over_sd sag_-100pA_mV sag_-100pA_over_sag_0pA sag_-100pA_over_2sd_0pA_in_step sag_ratio_0pA_own sag_0pA_over_hyper_-100pA sag_ratio_-100pA
1May2012_P3 0 0.814 2.393 8.364 17.539 7.330 10.776 3.146 0.061 0.448
20Sept2011_P2 0 0.592 1.348 2.848 7.936 5.889 6.706 0.707 0.038 0.224
24October2012_C2 0 1.027 1.916 20.136 8.666 4.522 4.220 4.943 0.046 0.208
Reading (b) on the 0 pA sweeps. sag_ratio_0pA_own divides the 0 pA 'sag' by the 0 pA sweep's own
deflection below baseline, which is noise too, so it is not a sag ratio of a hyperpolarization;
sag_0pA_over_hyper_-100pA divides it by the -100 pA hyperpolarization, the ratio noise alone would
add to the -100 pA sag ratio. Threshold of reading (b): 0.1
0 pA sweeps passing reading (b) on their own deflection: 3 of 3
noise ratio against the -100 pA hyperpolarization: [0.061, 0.038, 0.046] ; passing: 0 of 3
fig, axes = plt.subplots(1, 3, figsize=(13, 3.6), sharey=True)
for ax, (name, c) in zip(axes, cells.items()):
on, off = on_off[name]
for lvl, col in ((hc9.HYPER_STEP_PA, "tab:blue"), (hc9.DEPOL_STEP_PA, "tab:red")):
k = hc9.sweep_for(c, lvl)
ax.plot(c["t"] * 1e3, c["v"][:, k], color=col, lw=0.6, label=f"{lvl:+.0f} pA")
row = sag.set_index("cell").loc[name]
ax.axhline(row.v_min, color="0.5", ls=":", lw=0.8)
ax.axhline(row.v_end, color="0.5", ls="--", lw=0.8)
ax.set_title(f"{name}\n{panel[name]}", fontsize=9)
ax.set_xlabel("time (ms)")
axes[0].set_ylabel("membrane potential (mV)")
axes[0].legend(fontsize=8, loc="upper right")
plt.tight_layout()
plt.show()
yr = [(c["v"][:, hc9.sweep_for(c, l)].min(), c["v"][:, hc9.sweep_for(c, l)].max())
for c in cells.values() for l in (hc9.HYPER_STEP_PA, hc9.DEPOL_STEP_PA)]
print("Plotted: full 0 to 2000 ms sweeps; voltage range over the plotted traces",
round(min(a for a, _ in yr), 1), "to", round(max(b for _, b in yr), 1), "mV; nothing falls outside the axes.",
"Dotted line V_min, dashed line V_end (mean of the last 100 ms of the -100 pA step).")
Plotted: full 0 to 2000 ms sweeps; voltage range over the plotted traces -104.6 to 36.8 mV; nothing falls outside the axes. Dotted line V_min, dashed line V_end (mean of the last 100 ms of the -100 pA step).
The first 300 ms of the -100 pA step in each cell, the part of the trace where the sag appears, with the pre-step baseline, V_min and V_end marked, for comparison with the traces printed in Figure 7A and 7D of the paper (the figure is not reproduced here).
fig, axes = plt.subplots(1, 3, figsize=(13, 3.4), sharey=True)
ranges = []
for ax, (name, c) in zip(axes, cells.items()):
on, off = on_off[name]
k = hc9.sweep_for(c, hc9.HYPER_STEP_PA)
w = (c["t"] >= on - 0.05) & (c["t"] < on + 0.3)
x, y = (c["t"][w] - on) * 1e3, c["v"][w, k]
ranges.append((name, x.min(), x.max(), y.min(), y.max()))
ax.plot(x, y, color="tab:blue", lw=0.7)
row = sag.set_index("cell").loc[name]
ax.axhline(row.v_base, color="0.3", lw=0.8, label="baseline")
ax.axhline(row.v_min, color="0.5", ls=":", lw=0.8, label="V_min")
ax.axhline(row.v_end, color="0.5", ls="--", lw=0.8, label="V_end")
ax.axvline(row.t_min_ms, color="0.7", lw=0.6)
ax.set_title(f"{name}\n{panel[name]}, -100 pA", fontsize=9)
ax.set_xlabel("time from step onset (ms)")
axes[0].set_ylabel("membrane potential (mV)")
axes[0].legend(fontsize=8, loc="upper right")
plt.tight_layout()
plt.show()
for name, x0, x1, y0, y1 in ranges:
print(f"{name}: plotted {x0:.0f} to {x1:.0f} ms from onset, {y0:.1f} to {y1:.1f} mV; nothing outside the axes")
1May2012_P3: plotted -50 to 300 ms from onset, -102.8 to -62.9 mV; nothing outside the axes 20Sept2011_P2: plotted -50 to 300 ms from onset, -99.0 to -63.1 mV; nothing outside the axes 24October2012_C2: plotted -50 to 300 ms from onset, -104.6 to -62.8 mV; nothing outside the axes
print(f"Total run time of the notebook: {time.time() - started:.0f} s")
Total run time of the notebook: 7 s
Result: reproduced, resting on the measure the claims file assumes. On the -100 pA step all three cells show a sag under both readings, with V_min from a 1 ms running mean. (a) The relaxation from the peak hyperpolarization to the last 100 ms is 17.5, 7.9 and 8.7 mV for 1May2012_P3, 20Sept2011_P2 and 24October2012_C2, 22 to 47 times the threshold of twice the pre-step SD (0.19 to 0.56 mV). (b) The sag ratios are 0.448, 0.224 and 0.208, all above 0.1. Both readings give 3 of 3, and so does every variant shown (raw V_min, end windows of 50 to 200 ms): the smallest sag is 42.4 baseline SDs and the smallest ratio 0.200, so the outcome holds up to a 10 SD threshold or a ratio threshold of 0.15; at 0.2, two of three cells pass in every variant.
The 2 SD threshold of reading (a) does not separate a sag from noise in these recordings. On each cell's 0 pA sweep, with no current injected, the same measure gives 1.3 to 2.4 mV, 2.8 to 20.1 pre-step SDs, so all three 0 pA sweeps pass it; the 75 ms pre-step baseline misses the slow fluctuations of the voltage over the 1 s step (SD 0.59 to 1.03 mV in the 0 pA step window). On a synthetic trace without sag the raw minimum passes reading (a) in 100 of 100 draws (median 3.96 SD) and the 1 ms running mean in none (median 0.84 SD); reading (b) passes in none either way. Against the 0 pA sweep the -100 pA sags are 4.5 to 7.3 times larger, and 4.2 to 10.8 times twice the SD of the 0 pA step window. For reading (b), the 0 pA "sag" is 0.038 to 0.061 of each cell's -100 pA hyperpolarization, below the 0.1 threshold, while the -100 pA sag ratios are 0.208 to 0.448; reading (b) gives 3 of 3. Applied to the 0 pA sweep alone, reading (b) passes in all three cells (0.71 to 4.94), because with no hyperpolarization the ratio divides noise by noise; that number is not a false-positive level of the ratio on a hyperpolarizing step.
The cell the old description places in Figure 7D (24October2012_C2, P15, Sst-Cre;ChR2) has a sag ratio of 0.208, close to the 0.224 of 20Sept2011_P2 and below the 0.448 of 1May2012_P3, the two Figure 7A cells; that the 7D cell shows a sag too is the similarity the criterion judges. At +50 pA the three cells fire 15, 15 and 18 spikes in the 1 s step, which is shown only as a description. The paper defines no sag measure, so the outcome rests on the measure the claims file chose.
Summary¶
| Claim | Outcome |
|---|---|
| fig7a-olm-sag | reproduced (assumed measure of sag) |
| The other 25 claims (Figures 1 to 6 and 8 to 12, the light responses of Figure 7, the age analysis) | blocked, data missing |
How the analysis differs from the paper's. The paper shows the sag only in example traces
and defines no measure; the two readings of the claims file are used, and both agree. V_min is
taken from a 1 ms running mean, a standard choice. The 2 SD threshold of reading (a) is also
passed by each cell's 0 pA sweep, so that reading alone does not separate a sag from noise in
these recordings. On the 0 pA sweeps the noise "sag" is 0.038 to 0.061 of the -100 pA
hyperpolarization, below reading (b)'s 0.1, and the -100 pA sags are 4.5 to 7.3 times the 0 pA
value; on these references the outcome is the same. Which data
set cell is in which panel of Figure 7 comes from the superseded description in the bucket's
hc-9/old/ folder, not from the paper or the files; the criterion judges every cell alike, so
the outcome does not depend on it.
What was found about the data set:
- The data set holds only the current-step recordings of the three O-LM cells of Figure 7. None of the Cajal-Retzius recordings behind the paper's main findings are in it, nor the light-evoked responses of Figure 7C and 7D.
- The description gives the current steps as -100:30:+250 pA; the files hold 15 steps from -100 to +250 pA in 25 pA increments.
- The description says the injection usually starts about half a second into the sweep; in these files it starts at 75 ms and lasts 1 s.
- The current description does not say which animal, age or mouse line each cell comes from, or
which panel of Figure 7 shows it. The superseded description in
hc-9/old/does (P13 and P14 CXCR4-EGFP, P15 Sst-IRES-Cre;ChR2(H134R)-EYFP), and it matches the Figure 7 legends. That old description calls the files ABF; they are ATF text. - The zip's MD5 matches
checksums.md5; no recording is deposited twice; the voltage and current channels share one clock.