Gao, Ramachandran, Yu & He (2025) on dandiset 000637¶

Paper. Gao H, Ramachandran S, Yu K, He B (2025). Transcranial focused ultrasound modulates feedforward and feedback cortico-thalamo-cortical pathways by selectively activating excitatory neurons. J Neurosci 45(23):e2218242025, doi:10.1523/jneurosci.2218-24.2025.

Dataset. DANDI:000637, version 0.250829.1249 (doi:10.48324/dandi.000637/0.250829.1249), streamed. The paper cites no version of 000637 (its data statement names a different dandiset, 000946), so the most recent published version is pinned. The five published versions (0.250711.1756 to 0.250829.1249) are compared below: they hold the same files with the same digests, so the choice does not change the data.

Claims file. claims/000637/10.1523-jneurosci.2218-24.2025.yaml, written from the version of record (PMC12139594), which has no supplement and links no code.

What the dataset holds. 292 recordings from 25 rats, each with spike times of sorted units from one 32-channel probe in S1 and the trial start and stop times. The condition is coded only in the file's identifier. There are no spike waveforms, no POm, auditory cortex or CG1-targeted recordings, and of the duty-cycle series only the four conditions whose pulse duration is 100, 200 or 400 us, the three pulse durations the description lists.

Every claim of the paper is about units classified as RSUs or FSUs, and the paper classifies them from waveform shape (Methods, "Neuronal type classification"). Without waveforms none of the claims can be tested as the claims file states them, and most also need POm or AC units or conditions that are not here. All 31 claims are blocked. The notebook shows what the dataset does hold: the conditions, the checks of the metadata against the data, the trial clock, the counts of units under the paper's 2,000 spike rule, a check of the claims file's N2 reading on simulated units, and, as variants that decide no claim, the responses of all S1 units pooled across cell types.

Claim Outcome Main reason
unit-counts blocked data_missing: no waveforms for RSU/FSU classification; no POm units
fig4-s1-rsu-timelocked-3000, fig4-post-3000, fig4-30hz-s1-rsu-inhibition blocked data_missing: no cell types (and no POm for fig4-post-3000)
fig4-pom-rsu-timelocked-3000, fig4-30hz-pom-rsu-delayed, fig4-fsu-no-change blocked data_missing: no POm units, no cell types
fig5-sham-no-response blocked data_missing: no cell types, no POm; three unlabelled BASELINE recordings
fig5-cg1-control-no-response blocked data_missing: no CG1-targeted recordings
fig6a-s1-rsu-timelocked-high-prf, fig6a-s1-rsu-delayed, fig6b-s1-fsu blocked data_missing: no cell types
fig6c-pom-rsu-timelocked, fig6c-pom-rsu-delayed, fig6d-pom-fsu blocked data_missing: no POm units
fig7a to fig7d, fig8a, fig8b, fig9a, fig9b blocked data_missing: no cell types; most duty cycles absent; no POm for 7c and 7d
fig10a to fig10ef (dSTTC and STTC) blocked data_missing: no POm units
fig11b, fig11c, fig11d blocked data_missing: no auditory cortex units
In [1]:
import numpy as np, pandas as pd, requests, time
from scipy import stats
import matplotlib.pyplot as plt
import tfus
pd.set_option("display.width", 200); pd.set_option("display.max_columns", 40); pd.set_option("display.max_rows", 200)
print(tfus.DANDISET_ID, tfus.DANDISET_VERSION)
000637 0.250829.1249

Versions¶

The asset lists of the published versions, with each file's SHA-256 digest.

In [2]:
API = "https://api.dandiarchive.org/api/dandisets/000637/versions/"
versions = [v["version"] for v in requests.get(API, timeout=60).json()["results"] if v["version"] != "draft"]
def asset_list(v):
    out, url = {}, f"{API}{v}/assets/?page_size=1000&metadata=1"
    while url:
        r = requests.get(url, timeout=120).json()
        for a in r["results"]:
            out[a["path"]] = (a["metadata"].get("digest") or {}).get("dandi:sha2-256")
        url = r["next"]
    return out
lists = {v: asset_list(v) for v in versions}
ref = lists[tfus.DANDISET_VERSION]
for v, l in lists.items():
    print(v, len(l), "files;", "identical paths and digests to the pinned version" if l == ref else
          f"differs: {len(set(l.items()) ^ set(ref.items()))} path/digest pairs")
0.250711.1756 292 files; identical paths and digests to the pinned version
0.250711.1944 292 files; identical paths and digests to the pinned version
0.250828.1935 292 files; identical paths and digests to the pinned version
0.250829.1249 292 files; identical paths and digests to the pinned version

Loading¶

All 292 files are streamed; spike times and trial times are read into memory.

In [3]:
t0 = time.time()
sessions = tfus.load_all()
print(f"{len(sessions)} recordings loaded in {time.time() - t0:.0f} s")
meta = pd.DataFrame([{k: v for k, v in s.items() if k not in ("trial_start", "trial_stop", "units")} |
                     dict(n_units=len(s["units"]), n_trials=len(s["trial_start"]),
                          n_units_2000=sum(len(u) > tfus.MIN_SPIKES for u in s["units"]))
                     for s in sessions])
meta["dc_pct"] = tfus.duty_cycle(meta.prf, meta.pd_us)
292 recordings loaded in 307 s

Data check¶

In [4]:
print("Files per subject:"); print(meta.groupby("subject").size().to_string())
print("\nsubject_id field equals the identifier's subject:", bool((meta.subject_id == meta.subject).all()))
print("\nDistinct values across files:")
for col in ["description", "unit_columns", "trial_columns", "electrode_locations", "n_electrodes", "device", "acquisition", "processing"]:
    print(f"  {col}: {meta[col].astype(str).value_counts().to_dict()}")
print("\nRecordings by PRF (Hz) and pulse duration (us), stimulated recordings:")
print(pd.crosstab(meta.prf, meta.pd_us).to_string())
print("\nRecordings by input voltage:", meta.volt.value_counts(dropna=False).to_dict())
print("Ultrasound duration (ms):", meta.ud_ms.value_counts(dropna=False).to_dict())
print("Baseline recordings (no condition in the identifier):", meta.loc[meta.baseline, "identifier"].tolist())
print("Repeated recordings (_2, _3, _4):", meta.loc[meta.repeat > 1, "identifier"].tolist())
print("\nDuty cycles (PRF x PD) present:", sorted(meta.dc_pct.dropna().round(2).unique().tolist()))
print("\nTrials per recording:", meta.n_trials.value_counts().to_dict())
print("Units per recording: min", meta.n_units.min(), "median", meta.n_units.median(), "max", meta.n_units.max(),
      "; total", meta.n_units.sum(), "; with more than 2,000 spikes", meta.n_units_2000.sum())

iti = np.concatenate([np.diff(s["trial_start"]) for s in sessions])
dur = np.concatenate([s["trial_stop"] - s["trial_start"] for s in sessions])
print(f"\nIntertrial interval (s): min {iti.min():.3f}, 1st pct {np.percentile(iti, 1):.3f}, median {np.median(iti):.3f}, "
      f"99th pct {np.percentile(iti, 99):.3f}, max {iti.max():.3f}  (paper: 2.5 s with 10% jitter, 2.25 to 2.75 s)")
print(f"Trial duration stop - start (s): {np.unique(dur.round(4))}")
srt = all(np.all(np.diff(u) >= 0) for s in sessions for u in s["units"])
print("Spike times sorted in every unit:", srt)
allsp = np.concatenate([u for s in sessions for u in s["units"]])
print(f"Spike times (s): min {allsp.min():.4f}, max {allsp.max():.1f}")
last_trial = np.array([s["trial_stop"][-1] for s in sessions]); last_spike = np.array([max(u[-1] for u in s["units"] if len(u)) for s in sessions])
print(f"Last spike minus last trial stop (s): median {np.median(last_spike - last_trial):.1f}, max {np.max(last_spike - last_trial):.1f}")
print("Session start times:", meta.start_time.str[:10].value_counts().to_dict())

# Everything in the files, to look for waveforms, notes or stimulus records outside the fields read above
from collections import Counter
paths = Counter(p for s in sessions for p in s["h5_paths"] if not p.startswith("specifications/"))
print("\nHDF5 paths (outside specifications/) and the number of files that have each:")
for p_, n in sorted(paths.items()): print(f"  {p_}: {n}")
print("Units table description:", meta.units_description.value_counts().to_dict())
print("Text fields under general/ (key: value -> files):",
      Counter((k, v) for t in meta.general_text for k, v in t.items()))
print("Text fields under general/subject/, other than subject_id:",
      Counter((k, v) for t in meta.subject_text for k, v in t.items() if k != "subject_id"))
print("Electrode table columns:", meta.electrode_columns.astype(str).value_counts().to_dict())
print("Electrode z values:", meta.electrode_z.astype(str).value_counts().to_dict())
print("Electrode groups (name: location, description):", meta.electrode_groups.astype(str).value_counts().to_dict())
print("\nVoltage field of the identifier -> volts read:",
      meta[~meta.baseline].groupby("volt_field").volt.agg(["first", "size"]).to_dict("index"))
assert set(meta.loc[~meta.baseline, "volt_field"]) <= {"5", "7", "35", "55"}
Files per subject:
subject
BH401     8
BH444    11
BH445    12
BH446    14
BH448    10
BH452    11
BH453    13
BH454    13
BH455     5
BH456    11
BH457     6
BH458    18
BH459    13
BH460     7
BH461     6
BH462    13
BH463    13
BH464    17
BH465    13
BH466    13
BH467    13
BH468    12
BH469    14
BH470    13
BH477    13

subject_id field equals the identifier's subject: True

Distinct values across files:
  description: {'Neural Spikes for S1 tFUS': 292}
  unit_columns: {"['spike_times']": 292}
  trial_columns: {"['start_time', 'stop_time']": 292}
  electrode_locations: {"['unknown']": 292}
  n_electrodes: {'32': 292}
  device: {"['A1x32-Poly3-10mm-50-177-A32']": 292}
  acquisition: {'[]': 292}
  processing: {'[]': 292}

Recordings by PRF (Hz) and pulse duration (us), stimulated recordings:
pd_us   100.0  200.0  400.0  600.0
prf                               
30.0       18     26     17      0
300.0      19     26     18      0
1500.0     21     29     18      1
3000.0     20     30      0      0
4500.0     19     27      0      0

Recordings by input voltage: {5.0: 249, 5.5: 37, nan: 3, 7.0: 2, 3.5: 1}
Ultrasound duration (ms): {67.0: 289, nan: 3}
Baseline recordings (no condition in the identifier): ['BH458_BASELINE', 'BH460_BASELINE', 'BH461_BASELINE']
Repeated recordings (_2, _3, _4): ['BH457_1500_200_67_5_2', 'BH457_3000_200_67_5_2', 'BH458_3000_200_67_5_3', 'BH458_4500_200_67_5_2', 'BH458_3000_200_67_5_2', 'BH458_30_200_67_5_2', 'BH464_1500_200_67_5_3', 'BH464_4500_200_67_5_4', 'BH464_1500_200_67_5_2', 'BH464_300_200_67_5_3', 'BH464_30_200_67_5_2', 'BH464_4500_200_67_5_2', 'BH464_30_200_67_5_3', 'BH464_4500_200_67_5_3', 'BH464_300_200_67_5_2', 'BH469_1500_100_67_5_2']

Duty cycles (PRF x PD) present: [0.3, 0.6, 1.2, 3.0, 6.0, 12.0, 15.0, 30.0, 45.0, 60.0, 90.0]

Trials per recording: {500: 219, 230: 45, 425: 13, 375: 11, 328: 1, 400: 1, 344: 1, 252: 1}
Units per recording: min 7 median 36.0 max 69 ; total 11970 ; with more than 2,000 spikes 10798

Intertrial interval (s): min 2.250, 1st pct 2.255, median 2.501, 99th pct 2.745, max 2.750  (paper: 2.5 s with 10% jitter, 2.25 to 2.75 s)
Trial duration stop - start (s): [2.2]
Spike times sorted in every unit: True
Spike times (s): min 0.0007, max 1357.9
Last spike minus last trial stop (s): median 3.5, max 109.3
Session start times: {'2022-03-01': 292}

HDF5 paths (outside specifications/) and the number of files that have each:
  acquisition: 292
  analysis: 292
  file_create_date: 292
  general: 292
  general/devices: 292
  general/devices/array: 292
  general/extracellular_ephys: 292
  general/extracellular_ephys/electrodes: 292
  general/extracellular_ephys/electrodes/filtering: 292
  general/extracellular_ephys/electrodes/group: 292
  general/extracellular_ephys/electrodes/group_name: 292
  general/extracellular_ephys/electrodes/id: 292
  general/extracellular_ephys/electrodes/imp: 292
  general/extracellular_ephys/electrodes/location: 292
  general/extracellular_ephys/electrodes/x: 292
  general/extracellular_ephys/electrodes/y: 292
  general/extracellular_ephys/electrodes/z: 292
  general/extracellular_ephys/shank1: 292
  general/institution: 292
  general/intracellular_ephys: 292
  general/optogenetics: 292
  general/optophysiology: 292
  general/subject: 292
  general/subject/age: 292
  general/subject/description: 292
  general/subject/sex: 292
  general/subject/species: 292
  general/subject/subject_id: 292
  identifier: 292
  intervals: 292
  intervals/trials: 292
  intervals/trials/id: 292
  intervals/trials/start_time: 292
  intervals/trials/stop_time: 292
  processing: 292
  scratch: 292
  session_description: 292
  session_start_time: 292
  specifications: 292
  stimulus: 292
  stimulus/presentation: 292
  stimulus/templates: 292
  timestamps_reference_time: 292
  units: 292
  units/id: 292
  units/spike_times: 292
  units/spike_times_index: 292
Units table description: {'units table': 292}
Text fields under general/ (key: value -> files): Counter({('institution', 'Carnegie Mellon University'): 292})
Text fields under general/subject/, other than subject_id: Counter({('age', 'P24W'): 292, ('description', 'HSD:WI rat'): 292, ('sex', 'M'): 292, ('species', 'Rattus norvegicus'): 292})
Electrode table columns: {"['x', 'y', 'z', 'imp', 'location', 'filtering', 'group', 'group_name']": 292}
Electrode z values: {'[0.0]': 292}
Electrode groups (name: location, description): {"{'shank1': {'location': 'brain area', 'description': 'electrode group for shank1'}}": 292}

Voltage field of the identifier -> volts read: {'35': {'first': 3.5, 'size': 1}, '5': {'first': 5.0, 'size': 249}, '55': {'first': 5.5, 'size': 37}, '7': {'first': 7.0, 'size': 2}}

Checks of the description against the data, from the output above:

  • The description says "Each recording has 500 trials"; the counts per file are printed above, and those with fewer than 500 are a metadata error.
  • The description calls 100, 200 and 400 us "ultrasound durations (UDs)" and says "The pulse duration was kept at 67ms". The paper uses the opposite terms (UD 67 ms; PD 200 us in the PRF series). The identifiers carry 67 in the ultrasound-duration place for every stimulated recording, consistent with the paper's terms; the description has them swapped. One recording uses 600 us, which the description does not list.
  • The last field of the identifier is read as the input voltage (5 V in the paper's Methods); recordings at 3.5, 5.5 and 7 V are present and are not described.
  • The units carry no electrode, waveform or quality column, and all electrode locations are "unknown": the units cannot be assigned to S1 layers, to channels, or to cell types.
  • The intertrial intervals are within the paper's 2.25 to 2.75 s, consistent with trial start being the tFUS onset, which the next check tests.
  • No file holds anything beyond the fields listed in the tree above: acquisition, processing, analysis, scratch and stimulus are empty, the Units table has spike_times and its index only (no waveform_mean or other dataset outside colnames), general holds only the institution, the subject and the electrodes, with no experiment_description, notes, stimulus or protocol field, the one electrode group has location "brain area", and every electrode has z = 0. Nothing in any file names a target (S1, POm, CG1), a sham session or a cell type.
  • The voltage field takes only the values printed above (5, 7, 35, 55), so reading a two-digit field as tenths of a volt affects only 35 and 55 (3.5 and 5.5 V).

Trial clock¶

The files do not say that a trial's start is the tFUS onset. For each unit with more than 2,000 spikes in a stimulated recording at 1500, 3000 or 4500 Hz (where the paper reports responses), the spike counts in the 67 ms windows after the trial start are compared with the 67 ms before it, at the true trial starts and at starts shifted by 1.2 s (half an intersonication interval). A response that appears at the true starts and not at the shifted ones shows the trials and spikes share a clock and that the start is the onset to within the window.

In [5]:
units = tfus.unit_table(sessions)
shifted = tfus.unit_table(sessions, shift_s=tfus.SHIFT_S)
hi = units.prf.isin([1500, 3000, 4500]) & ~units.baseline
hs = shifted.prf.isin([1500, 3000, 4500]) & ~shifted.baseline
rows = []
for name, tab, m in [("true onsets", units, hi), ("shifted 1.2 s", shifted, hs)]:
    for w in ["0-67", "67-134", "134-201", "201-268"]:
        r = stats.wilcoxon(tab.loc[m, f"count_{w}"], tab.loc[m, "count_pre"])
        rows.append(dict(onsets=name, window=w, n_units=int(m.sum()),
                         median_N1=tab.loc[m, f"n1_{w}"].median(),
                         frac_units_up=np.mean(tab.loc[m, f"count_{w}"] > tab.loc[m, "count_pre"]),
                         signed_rank_p=r.pvalue))
clock = pd.DataFrame(rows)
with pd.option_context("display.float_format", "{:.3g}".format):
    print(clock.to_string(index=False))
# Per subject: median N1 at 201-268 ms with true and shifted onsets
per = pd.DataFrame({"true": units[hi].groupby("subject")["n1_201-268"].median(),
                    "shifted": shifted[hs].groupby("subject")["n1_201-268"].median()})
print("\nPer subject, median N1 at 201-268 ms (%), PRF 1500-4500 Hz:"); print(per.round(1).T.to_string())
print("Subjects with true > shifted:", int((per["true"] > per["shifted"]).sum()), "of", len(per))
       onsets  window  n_units  median_N1  frac_units_up  signed_rank_p
  true onsets    0-67     6093        111          0.744              0
  true onsets  67-134     6093        102          0.578       1.21e-67
  true onsets 134-201     6093        106          0.689              0
  true onsets 201-268     6093        113          0.781              0
shifted 1.2 s    0-67     6093       98.6          0.431       6.21e-31
shifted 1.2 s  67-134     6093       97.9          0.398       9.08e-65
shifted 1.2 s 134-201     6093       97.7          0.392        6.2e-68
shifted 1.2 s 201-268     6093       97.6           0.39       1.95e-76

Per subject, median N1 at 201-268 ms (%), PRF 1500-4500 Hz:
subject  BH401  BH444  BH445  BH446  BH448  BH452  BH453  BH454  BH455  BH456  BH457  BH458  BH459  BH460  BH461  BH462  BH463  BH464  BH465  BH466  BH467  BH468  BH469  BH470  BH477
true     100.0   98.3  104.6  125.5  104.1  158.8  106.3   99.6  106.1  131.6  134.3  114.8  134.7   99.0  116.8  187.7  311.4  102.9  113.5  104.1  105.5  107.6  107.5  111.2  110.5
shifted   98.9  101.3  103.5   96.5  102.4   92.8   98.6  100.0  102.0   98.7   92.5   97.7   96.5   99.0  100.5   87.9   88.5   99.5   97.6   99.8   98.4   97.8   98.6   96.0   98.8
Subjects with true > shifted: 23 of 25
In [6]:
# Time course over the whole trial: per unit, rate in 67 ms bins from 1.2 s before to 2.28 s after the trial
# start, as a percentage of the unit's mean over all bins, then averaged over units
edges = np.arange(-18, 35) * tfus.WIN_MS / 1000.0
def psth_table(sess):
    r = tfus.population_psth(sess, edges)
    r = r[r.mean(1) > 0]
    return 100 * r / r.mean(1, keepdims=True), len(r)
hi_sess = [s for s in sessions if not s["baseline"] and s["prf"] in (1500, 3000, 4500)]
bl_sess = [s for s in sessions if s["baseline"]]
(ph, nh), (pb, nb_) = psth_table(hi_sess), psth_table(bl_sess)
psth = pd.DataFrame({"bin_start_ms": np.round(edges[:-1] * 1000).astype(int),
                     f"PRF 1500-4500 Hz (n={nh})": ph.mean(0).round(1),
                     f"BASELINE (n={nb_})": pb.mean(0).round(1)})
print(psth.to_string(index=False))
fig, ax = plt.subplots(figsize=(10, 3))
for col in psth.columns[1:]:
    ax.step(psth.bin_start_ms, psth[col], where="post", label=col)
ax.axvline(0, color="grey", lw=0.8); ax.axvline(1200, color="grey", lw=0.8, ls=":")
ax.set_xlabel("ms from trial start (bin start)"); ax.set_ylabel("% of unit mean"); ax.legend(); plt.tight_layout(); plt.show()
 bin_start_ms  PRF 1500-4500 Hz (n=6093)  BASELINE (n=119)
        -1206                       94.4              98.2
        -1139                       93.0              98.8
        -1072                       92.9              97.5
        -1005                       92.6              90.9
         -938                       92.8              90.4
         -871                       91.5             100.3
         -804                       91.0             100.2
         -737                       90.2              99.9
         -670                       90.5             101.5
         -603                       91.2             103.2
         -536                       90.4              95.3
         -469                       91.7              96.3
         -402                       91.2             100.8
         -335                       92.1             108.4
         -268                       93.1             104.7
         -201                       91.9             105.5
         -134                       92.8             101.5
          -67                       93.2              99.1
            0                      142.8              96.9
           67                       95.7             101.3
          134                      108.6             104.8
          201                      147.3             103.9
          268                      139.0             105.7
          335                      128.5              98.3
          402                      125.7             102.2
          469                      121.1             106.9
          536                      117.8             106.5
          603                      115.0             102.2
          670                      112.2              97.7
          737                      109.2              95.9
          804                      106.9              96.7
          871                      104.9              98.8
          938                      101.8              98.7
         1005                       99.0              97.6
         1072                       98.2              96.6
         1139                       96.5             102.5
         1206                       94.5              99.6
         1273                       93.6             100.0
         1340                       93.6              99.5
         1407                       93.0              96.2
         1474                       92.8             100.5
         1541                       91.8              95.3
         1608                       91.6             100.8
         1675                       91.3              96.3
         1742                       91.4             101.3
         1809                       91.3              99.6
         1876                       91.6              98.7
         1943                       90.8             103.2
         2010                       91.3              97.5
         2077                       91.4             101.9
         2144                       91.4             102.9
         2211                       92.1             101.0
No description has been provided for this image

Result: trial clock. At the true trial starts the response is a clear contrast with the shifted starts (table above), so the trial start is taken as the tFUS onset. The shifted onsets are not a null, though: with onsets moved by 1.2 s every window is below the shifted pre-window. The time course over the whole trial (PSTH above) shows why: after the onset the rate rises in 0 to 67 ms (143% of the unit mean), dips at 67 to 134 ms, peaks again at 201 to 268 ms (147%) and then decays slowly, reaching the unit mean only about 1 s after the onset and settling at 91 to 94% from about 1.3 s to the next trial. The shifted pre-window (1,133 to 1,200 ms) falls on the end of that decay (about 96%), and the shifted windows (1,200 to 1,468 ms) a little later (93 to 95%), so the shifted windows sit below their pre-window because the post-stimulus elevation lasts more than a second. This is a slow after-effect of the stimulus, not a sign that trial start and onset are misaligned: the rise at 0 to 67 ms is confined to the first bin after the trial start. The 1.2 s shift is therefore a contrast, not a null with no stimulus effect; the conclusion that trial start is the tFUS onset rests on the contrast, which holds in every window (pooled over units) and, at 201 to 268 ms, the one window compared per rat, in 23 of 25 rats.

Recordings deposited more than once¶

Each unit is fingerprinted by its spike count and its first 20 spike times (to 0.1 ms). A fingerprint shared by units in two files, or two units of one file, would be a recording deposited twice.

In [7]:
fp = {}
for s in sessions:
    for i, u in enumerate(s["units"]):
        fp.setdefault((len(u), tuple(np.round(u[:20], 4))), []).append((s["identifier"], i))
shared = [v for v in fp.values() if len(v) > 1]
print("Units:", sum(len(v) for v in fp.values()), "; fingerprints shared by more than one unit:", len(shared))
for v in shared[:10]: print("  ", v)
# Recordings with the same condition in the same rat (_2, _3): their unit counts
for subj, g in meta[meta.repeat > 1].groupby("subject"):
    print(subj, "repeats:", g.identifier.tolist(), "units", g.n_units.tolist())
Units: 11970 ; fingerprints shared by more than one unit: 0
BH457 repeats: ['BH457_1500_200_67_5_2', 'BH457_3000_200_67_5_2'] units [32, 32]
BH458 repeats: ['BH458_3000_200_67_5_3', 'BH458_4500_200_67_5_2', 'BH458_3000_200_67_5_2', 'BH458_30_200_67_5_2'] units [63, 60, 65, 61]
BH464 repeats: ['BH464_1500_200_67_5_3', 'BH464_4500_200_67_5_4', 'BH464_1500_200_67_5_2', 'BH464_300_200_67_5_3', 'BH464_30_200_67_5_2', 'BH464_4500_200_67_5_2', 'BH464_30_200_67_5_3', 'BH464_4500_200_67_5_3', 'BH464_300_200_67_5_2'] units [34, 35, 32, 32, 32, 32, 32, 31, 31]
BH469 repeats: ['BH469_1500_100_67_5_2'] units [32]

What the claims need, and what is here¶

The claims test cell-typed units (RSU, FSU) in S1 and POm, and for Figure 11 in auditory cortex, under the paper's conditions. The table lists each condition set the claims use, with the pulse durations that give it (PD = DC / PRF, Figure 1B legend), and how many recordings of each the dataset has.

In [8]:
needed = {
    "PRF series, 200 us (Fig 4, 6, 10C)": [(p, 200) for p in (30, 300, 1500, 3000, 4500)],
    "3000 Hz DC series 0.6-90% (Fig 7, 10F)": [(3000, d) for d in (2, 20, 100, 200, 300)],
    "30 Hz DC series 0.6-90% (Fig 8, 10E)": [(30, d) for d in (200, 2000, 10000, 20000, 30000)],
    "60% DC PRF series (Fig 9, 10D)": [(30, 20000), (300, 2000), (1500, 400), (3000, 200), (4500, 133)],
}
stim = meta[~meta.baseline]
rows = []
for name, conds in needed.items():
    for prf, pd_us in conds:
        n = int(((stim.prf == prf) & (stim.pd_us == pd_us)).sum())
        rows.append(dict(series=name, prf=prf, pd_us=pd_us, dc_pct=round(prf * pd_us * 1e-4, 2), recordings=n))
print(pd.DataFrame(rows).to_string(index=False))
print("\nSham (no tFUS): no recording is labelled sham; recordings named BASELINE:", int(meta.baseline.sum()))
print("CG1-targeted control: no identifier or field names a target; S1 is the only target described.")
print("Brain areas: one probe per file, electrode locations", meta.electrode_locations.astype(str).unique().tolist(),
      "; the description names only S1. No POm or auditory cortex units.")
print("Waveforms or waveform features in the Units table:", any(c != ["spike_times"] for c in meta.unit_columns))
                                series  prf  pd_us  dc_pct  recordings
    PRF series, 200 us (Fig 4, 6, 10C)   30    200    0.60          26
    PRF series, 200 us (Fig 4, 6, 10C)  300    200    6.00          26
    PRF series, 200 us (Fig 4, 6, 10C) 1500    200   30.00          29
    PRF series, 200 us (Fig 4, 6, 10C) 3000    200   60.00          30
    PRF series, 200 us (Fig 4, 6, 10C) 4500    200   90.00          27
3000 Hz DC series 0.6-90% (Fig 7, 10F) 3000      2    0.60           0
3000 Hz DC series 0.6-90% (Fig 7, 10F) 3000     20    6.00           0
3000 Hz DC series 0.6-90% (Fig 7, 10F) 3000    100   30.00          20
3000 Hz DC series 0.6-90% (Fig 7, 10F) 3000    200   60.00          30
3000 Hz DC series 0.6-90% (Fig 7, 10F) 3000    300   90.00           0
  30 Hz DC series 0.6-90% (Fig 8, 10E)   30    200    0.60          26
  30 Hz DC series 0.6-90% (Fig 8, 10E)   30   2000    6.00           0
  30 Hz DC series 0.6-90% (Fig 8, 10E)   30  10000   30.00           0
  30 Hz DC series 0.6-90% (Fig 8, 10E)   30  20000   60.00           0
  30 Hz DC series 0.6-90% (Fig 8, 10E)   30  30000   90.00           0
        60% DC PRF series (Fig 9, 10D)   30  20000   60.00           0
        60% DC PRF series (Fig 9, 10D)  300   2000   60.00           0
        60% DC PRF series (Fig 9, 10D) 1500    400   60.00          18
        60% DC PRF series (Fig 9, 10D) 3000    200   60.00          30
        60% DC PRF series (Fig 9, 10D) 4500    133   59.85           0

Sham (no tFUS): no recording is labelled sham; recordings named BASELINE: 3
CG1-targeted control: no identifier or field names a target; S1 is the only target described.
Brain areas: one probe per file, electrode locations ["['unknown']"] ; the description names only S1. No POm or auditory cortex units.
Waveforms or waveform features in the Units table: False

unit-counts (blocked): counts under the 2,000 spike rule¶

The analyses of PRF at constant 200 us pulse duration rest on 1,050 S1 RSUs, 247 S1 FSUs, 913 POm RSUs and 106 POm FSUs, and the 3000 Hz and 30 Hz sessions of Figure 4 on 312 and 148 S1 RSUs, 54 and 61 S1 FSUs, 233 and 189 POm RSUs, and 17 and 43 POm FSUs.

The RSU/FSU split needs waveforms, which the dataset lacks, and the POm counts need POm units. As a variant that does not decide the claim, the S1 units with more than 2,000 spikes are counted per condition, all cell types together, and compared with the paper's S1 RSU + FSU sums. The counts are given for all input voltages and for 5 V alone (the paper's stated input voltage).

In [9]:
s1 = units[~units.baseline]
def count(df, prf=None, pd_us=200):
    m = (df.pd_us == pd_us) & (df.prf.isin(prf) if prf is not None else True)
    return int(m.sum())
paper = {"PRF series, all PRFs": 1050 + 247, "3000 Hz 200 us": 312 + 54, "30 Hz 200 us": 148 + 61}
rows = []
for label, prfs in [("PRF series, all PRFs", None), ("3000 Hz 200 us", [3000]), ("30 Hz 200 us", [30])]:
    rows.append(dict(set=label, paper_S1_RSU_plus_FSU=paper[label],
                     dataset_all_volts=count(s1, prfs), dataset_5V=count(s1[s1.volt == 5], prfs)))
print(pd.DataFrame(rows).to_string(index=False))
print("\nS1 units with >2,000 spikes per PRF at 200 us (all voltages):", s1[s1.pd_us == 200].groupby("prf").size().to_dict())
print("Paper's S1 RSUs per PRF at 200 us (Results, Fig 6A): 1500 Hz 284, 3000 Hz 312, 4500 Hz 152; S1 FSUs at 60% DC 1500 Hz 48 (Fig 9B)")
print(f"Ratio of dataset to paper, PRF series all voltages: {rows[0]['dataset_all_volts'] / rows[0]['paper_S1_RSU_plus_FSU']:.2f}")
print("Units dropped by the 2,000 spike rule:", int(meta.n_units.sum() - meta.n_units_2000.sum()), "of", int(meta.n_units.sum()))
# Subjects of 000637 and of 000946, the dandiset the paper's data statement names, at their pinned versions
def subjects(ds, v):
    r = requests.get(f"https://api.dandiarchive.org/api/dandisets/{ds}/versions/{v}/assets/paths/?path_prefix=&page_size=1000",
                     timeout=60).json()
    assert r["next"] is None
    return sorted(x["path"] for x in r["results"])
s637, s946 = subjects("000637", tfus.DANDISET_VERSION), subjects("000946", "0.250829.1249")
print("\n000637 top-level folders:", len(s637), s637)
print("000946 (version 0.250829.1249) top-level folders:", len(s946), s946)
rats946 = {x.replace("sub-STTC-", "sub-") for x in s946}
print("000946 rats:", len(rats946), min(rats946), "to", max(rats946), "; shared with 000637:", sorted(rats946 & set(s637)))
                 set  paper_S1_RSU_plus_FSU  dataset_all_volts  dataset_5V
PRF series, all PRFs                   1297               4642        4023
      3000 Hz 200 us                    366               1088         902
        30 Hz 200 us                    209                862         767

S1 units with >2,000 spikes per PRF at 200 us (all voltages): {30.0: 862, 300.0: 853, 1500.0: 938, 3000.0: 1088, 4500.0: 901}
Paper's S1 RSUs per PRF at 200 us (Results, Fig 6A): 1500 Hz 284, 3000 Hz 312, 4500 Hz 152; S1 FSUs at 60% DC 1500 Hz 48 (Fig 9B)
Ratio of dataset to paper, PRF series all voltages: 3.58
Units dropped by the 2,000 spike rule: 1172 of 11970
000637 top-level folders: 25 ['sub-BH401', 'sub-BH444', 'sub-BH445', 'sub-BH446', 'sub-BH448', 'sub-BH452', 'sub-BH453', 'sub-BH454', 'sub-BH455', 'sub-BH456', 'sub-BH457', 'sub-BH458', 'sub-BH459', 'sub-BH460', 'sub-BH461', 'sub-BH462', 'sub-BH463', 'sub-BH464', 'sub-BH465', 'sub-BH466', 'sub-BH467', 'sub-BH468', 'sub-BH469', 'sub-BH470', 'sub-BH477']
000946 (version 0.250829.1249) top-level folders: 39 ['sub-BH494', 'sub-BH495', 'sub-BH499', 'sub-BH500', 'sub-BH501', 'sub-BH510', 'sub-BH511', 'sub-BH513', 'sub-BH528', 'sub-BH530', 'sub-BH532', 'sub-BH535', 'sub-BH537', 'sub-BH538', 'sub-BH539', 'sub-BH547', 'sub-BH551', 'sub-BH553', 'sub-BH558', 'sub-BH594', 'sub-STTC-BH494', 'sub-STTC-BH495', 'sub-STTC-BH499', 'sub-STTC-BH500', 'sub-STTC-BH501', 'sub-STTC-BH510', 'sub-STTC-BH511', 'sub-STTC-BH513', 'sub-STTC-BH528', 'sub-STTC-BH530', 'sub-STTC-BH532', 'sub-STTC-BH535', 'sub-STTC-BH537', 'sub-STTC-BH538', 'sub-STTC-BH539', 'sub-STTC-BH547', 'sub-STTC-BH551', 'sub-STTC-BH553', 'sub-STTC-BH558']
000946 rats: 20 sub-BH494 to sub-BH594 ; shared with 000637: []

Result: unit-counts, blocked. The six POm counts cannot be checked (no POm units) and the six S1 counts need the RSU/FSU split. Even so the S1 totals rule out "reproduced": whatever classification were applied, each S1 RSU + FSU pair would add up to the dataset's total, 4,642, 1,088 and 862 units against the paper's 1,297, 366 and 209, so at most one count of each S1 pair could match, and at most 6 + 3 = 9 of the 12 counts, below the criterion's 10. The outcome is kept as blocked because the POm half cannot be computed.

Whether 000637 holds the recordings the paper analysed at all is open, and the numbers point against it being only those. The paper's data statement names DANDI:000946, not 000637, and 000946 holds 20 other rats (BH494 to BH594, none shared with 000637's BH401 to BH477; folder lists printed above). The paper reports two groups of 16 and 12 rats (28), and 000637 has 25. The dataset has 3.6 times as many S1 units over the paper's 2,000 spike threshold in the PRF series as the paper's S1 RSUs and FSUs together (4,023 at 5 V alone, 3.1 times). 000637 may hold more sessions or rats than the paper used, the S1 probes of other groups, or different recordings; nothing in the files says which (data check).

Variant: all S1 units pooled, PRF series at 200 us¶

This decides no claim: the claims are about RSUs and FSUs separately, and the paper's RSUs and FSUs respond differently. It shows whether the S1 recordings in this dataset carry the time-locked and delayed responses that the paper reports for S1 units. Normalized rates N1 and N2 are those of the claims file's Method S; the tests are Wilcoxon signed-rank (W2) and rank-sum (W1) of the raw rates in each window against the pre-stimulus window, one-tailed for an increase and for a decrease, Bonferroni over the four windows (M1). Units are counted once per recording.

In [10]:
def window_tests(df, label):
    rows = []
    for w in ["0-67", "67-134", "134-201", "201-268"]:
        x, b = df[f"rate_{w}"], df["rate_pre"]
        r = dict(set=label, window=w, n=len(df), mean_N1=df[f"n1_{w}"].mean(), sem_N1=df[f"n1_{w}"].sem(),
                 mean_N2=df[f"n2_{w}"].mean())
        for tail, alt in [("up", "greater"), ("down", "less")]:
            r[f"W2_{tail}_p"] = min(1.0, 4 * stats.wilcoxon(x, b, alternative=alt).pvalue)
            r[f"W1_{tail}_p"] = min(1.0, 4 * stats.mannwhitneyu(x, b, alternative=alt).pvalue)
        rows.append(r)
    return rows
prf_rows = []
for prf in (30, 300, 1500, 3000, 4500):
    prf_rows += window_tests(s1[(s1.pd_us == 200) & (s1.prf == prf)], f"{prf} Hz")
prf_tab = pd.DataFrame(prf_rows)
with pd.option_context("display.float_format", "{:.3g}".format):
    print(prf_tab.to_string(index=False))
for w in ["0-67", "201-268"]:
    groups = [s1.loc[(s1.pd_us == 200) & (s1.prf == p), f"n1_{w}"].dropna() for p in (30, 300, 1500, 3000, 4500)]
    print(f"Kruskal-Wallis across PRFs, N1 {w}: H = {stats.kruskal(*groups).statistic:.1f}, p = {stats.kruskal(*groups).pvalue:.3g}")
print("Units with no pre-stimulus spike (N1 undefined):", int(s1['n1_0-67'].isna().sum()))
print("Range of N1 at 201-268 ms (%):", round(s1['n1_201-268'].min(), 1), "to", round(s1['n1_201-268'].max(), 1))
    set  window    n  mean_N1  sem_N1  mean_N2   W2_up_p  W1_up_p  W2_down_p  W1_down_p
  30 Hz    0-67  862      105   0.428     69.2  2.57e-26  0.00667          1          1
  30 Hz  67-134  862      101   0.303     65.1     0.548        1          1          1
  30 Hz 134-201  862     99.2   0.305     64.5         1        1     0.0785      0.718
  30 Hz 201-268  862     98.2   0.336     63.6         1        1   3.78e-05      0.275
 300 Hz    0-67  853      126    2.11     83.3  6.84e-65  1.1e-18          1          1
 300 Hz  67-134  853      100   0.534     63.9         1        1    0.00555          1
 300 Hz 134-201  853      105   0.984     66.3     0.989   0.0953          1          1
 300 Hz 201-268  853      132    4.28     69.8  3.44e-09 0.000334          1          1
1500 Hz    0-67  938      231    29.5      108  4.37e-96 1.21e-45          1          1
1500 Hz  67-134  938      101   0.363     63.2  0.000665    0.931          1          1
1500 Hz 134-201  938      121    2.67     65.8  3.16e-49 5.85e-05          1          1
1500 Hz 201-268  938      231    21.7     72.7  2.54e-52 2.24e-13          1          1
3000 Hz    0-67 1088      364     189     85.5  2.69e-87 7.17e-25          1          1
3000 Hz  67-134 1088      101   0.541     60.5     0.435        1          1          1
3000 Hz 134-201 1088      112    1.09     63.2  2.17e-39 0.000739          1          1
3000 Hz 201-268 1088      132    1.96     68.4 9.38e-107 3.19e-17          1          1
4500 Hz    0-67  901      106   0.667     67.1  8.48e-28  0.00155          1          1
4500 Hz  67-134  901      109   0.726     65.8  6.16e-30  0.00221          1          1
4500 Hz 134-201  901      140    3.89     69.5  3.96e-68 3.34e-14          1          1
4500 Hz 201-268  901      241    25.7     78.9  1.29e-78 5.67e-32          1          1
Kruskal-Wallis across PRFs, N1 0-67: H = 318.2, p = 1.27e-67
Kruskal-Wallis across PRFs, N1 201-268: H = 660.9, p = 1.01e-141
Units with no pre-stimulus spike (N1 undefined): 31
Range of N1 at 201-268 ms (%): 38.5 to 89700.0
In [11]:
# By animal: median N1 per PRF at 0-67 and 201-268 ms (S1 units, 200 us)
d = s1[s1.pd_us == 200]
for w in ["0-67", "201-268"]:
    t = d.pivot_table(index="subject", columns="prf", values=f"n1_{w}", aggfunc="median")
    print(f"Median N1 {w} ms (%), by subject and PRF:"); print(t.round(0).to_string())
    print("Subjects with median above 100% at each PRF:", (t > 100).sum().to_dict(), "of", t.notna().sum().to_dict(), "\n")
fig, ax = plt.subplots(1, 2, figsize=(10, 3.5), sharey=True)
for a, w in zip(ax, ["0-67", "201-268"]):
    sub = prf_tab[prf_tab.window == w]
    a.errorbar(range(5), sub.mean_N1, yerr=sub.sem_N1, fmt="o-", color="k")
    a.axhline(100, color="grey", lw=0.8); a.set_xticks(range(5)); a.set_xticklabels(sub.set)
    a.set_title(f"All S1 units, {w} ms, PD 200 us"); a.set_ylabel("mean N1 (%) +/- SEM")
plt.tight_layout(); plt.show()
print("Plotted means (%):", prf_tab[prf_tab.window.isin(['0-67', '201-268'])].mean_N1.round(1).tolist())
Median N1 0-67 ms (%), by subject and PRF:
prf      30.0    300.0   1500.0  3000.0  4500.0
subject                                        
BH401       NaN   112.0   151.0   104.0     NaN
BH444     104.0    96.0   105.0    99.0    95.0
BH445      92.0   133.0   258.0   198.0   104.0
BH446     110.0    92.0   104.0   101.0    87.0
BH448      99.0   100.0   110.0   103.0    94.0
BH452     100.0    99.0    97.0   115.0   111.0
BH453     101.0   102.0   112.0   100.0   102.0
BH454     100.0   115.0   130.0   101.0   105.0
BH455     101.0   108.0   114.0   192.0   104.0
BH456     102.0   103.0   139.0   117.0     NaN
BH457       NaN     NaN    99.0   101.0    98.0
BH458     113.0   120.0   110.0   108.0   107.0
BH459     105.0   103.0   105.0   134.0   103.0
BH460     106.0   116.0   144.0   149.0   107.0
BH461     106.0    96.0   117.0   107.0   103.0
BH462     106.0   107.0   148.0   107.0   104.0
BH463     102.0   136.0   103.0   123.0   109.0
BH464     104.0   109.0   126.0   111.0   103.0
BH465     100.0   113.0   133.0   123.0   104.0
BH466      99.0   109.0   122.0   130.0   103.0
BH467     110.0   126.0   187.0   164.0   110.0
BH468     115.0   173.0   106.0   175.0   122.0
BH469      99.0   109.0   128.0   114.0   104.0
BH470     101.0   145.0   338.0   212.0   111.0
BH477     108.0   130.0   303.0   125.0   109.0
Subjects with median above 100% at each PRF: {30.0: 17, 300.0: 19, 1500.0: 23, 3000.0: 24, 4500.0: 19} of {30.0: 23, 300.0: 24, 1500.0: 25, 3000.0: 25, 4500.0: 23} 

Median N1 201-268 ms (%), by subject and PRF:
prf      30.0    300.0   1500.0  3000.0  4500.0
subject                                        
BH401       NaN    91.0   100.0    98.0     NaN
BH444      92.0   109.0    94.0   104.0   100.0
BH445      76.0    83.0    97.0   106.0   103.0
BH446     104.0   335.0   108.0    99.0   519.0
BH448      98.0    98.0   101.0   100.0    96.0
BH452      86.0   100.0   127.0   162.0   255.0
BH453      98.0    99.0   103.0    98.0   108.0
BH454     104.0   101.0   100.0    98.0    99.0
BH455     101.0   100.0    98.0   139.0   104.0
BH456      98.0   101.0   117.0   104.0     NaN
BH457       NaN     NaN   143.0   125.0   157.0
BH458      96.0   104.0   151.0   115.0   116.0
BH459     103.0   104.0   114.0   146.0   182.0
BH460      95.0   105.0   103.0    90.0    95.0
BH461     101.0    96.0   113.0   164.0   108.0
BH462     104.0   101.0   140.0   120.0   181.0
BH463     101.0   187.0   444.0   119.0   256.0
BH464     101.0   100.0   101.0   104.0   104.0
BH465      96.0    99.0   103.0   130.0   112.0
BH466     102.0   105.0   100.0   101.0   109.0
BH467     104.0   106.0   102.0    98.0    98.0
BH468     101.0   100.0    98.0   112.0   105.0
BH469      98.0   103.0   104.0   103.0   101.0
BH470     101.0    99.0   102.0   100.0   154.0
BH477     100.0    99.0   109.0   104.0   100.0
Subjects with median above 100% at each PRF: {30.0: 12, 300.0: 13, 1500.0: 20, 3000.0: 18, 4500.0: 18} of {30.0: 23, 300.0: 24, 1500.0: 25, 3000.0: 25, 4500.0: 23} 

No description has been provided for this image
Plotted means (%): [104.9, 98.2, 125.8, 131.8, 230.6, 230.5, 364.4, 132.3, 106.2, 241.1]

Check of the two normalized-rate readings on units with no response¶

Method S leaves open N1 (summed counts over summed pre-stimulus counts) and N2 (mean over trials of the per-trial ratio, trials without a pre-stimulus spike left out). Poisson units with no response, 500 trials, show what each reading gives when nothing happens.

In [12]:
rng = np.random.default_rng(0)
fake = []
for rate in (3, 10, 30, 100):
    onsets = np.arange(500) * 2.5 + 1.0
    us = [np.sort(rng.uniform(0, 1252.0, rng.poisson(rate * 1252.0))) for _ in range(40)]
    fake.append(dict(identifier=f"poisson_{rate}Hz", subject="sim", prf=rate, pd_us=0, volt=0, baseline=False,
                     trial_start=onsets, units=us))
ft = tfus.unit_table(fake)
cols = [f"n{k}_{w}" for k in (1, 2) for w in ["0-67", "201-268"]]
print("Mean normalized rate (%) of 40 Poisson units with no response, by firing rate (Hz):")
print(ft.groupby("prf")[cols].mean().round(1).to_string())
print("\nMean N2 (%) of the dataset's S1 units, PRF series at 200 us, and of BASELINE units, with their median rate (Hz):")
for lab, d_ in [("PRF series", s1[s1.pd_us == 200]), ("BASELINE", units[units.baseline])]:
    print(f"  {lab}: N2 0-67 {d_['n2_0-67'].mean():.1f}, N2 201-268 {d_['n2_201-268'].mean():.1f}, "
          f"median pre-stimulus rate {d_['rate_pre'].median():.1f}")
Mean normalized rate (%) of 40 Poisson units with no response, by firing rate (Hz):
     n1_0-67  n1_201-268  n2_0-67  n2_201-268
prf                                          
3      101.3       101.2     18.8        18.7
10     100.2       100.4     56.0        57.8
30     100.5       100.3    115.6       115.7
100     99.6        99.7    121.2       121.2

Mean N2 (%) of the dataset's S1 units, PRF series at 200 us, and of BASELINE units, with their median rate (Hz):
  PRF series: N2 0-67 83.1, N2 201-268 70.7, median pre-stimulus rate 12.0
  BASELINE: N2 0-67 63.7, N2 201-268 59.3, median pre-stimulus rate 11.5

Result: N2 reading. N1 is unbiased on Poisson units with no response (99.6 to 101.3% at every rate). N2 is strongly biased, and the sign of the bias depends on the firing rate: 18.8% at 3 Hz, 56 to 58% at 10 Hz, 116% at 30 Hz and 121% at 100 Hz. A per-trial ratio with a pre-stimulus count of 0 to a few spikes in 67 ms is dominated by trials where the window count is 0 (pulling the mean down) or where the pre-stimulus count is 1 (pulling it up). The dataset's S1 units, whose median pre-stimulus rate is 12 Hz, have mean N2 between 60 and 86% in every window but one (108% at 1500 Hz, 0 to 67 ms, where mean N1 is 231%), including the BASELINE recordings (59 to 64%) and windows where N1 is above 200%. Under reading R2 (normalized rates tested against 100%), N2 would report decreases in units with no response and mask increases, so a claim tested on it could be blocked or not reproduced for an artefactual reason. No claim is tested here, so this changes no outcome; it is asked of a reader in questions/000637.yaml.

Variant: BASELINE recordings¶

Three recordings are named BASELINE, with 500 trials each (printed below). Nothing in the files says whether they are the paper's sham sessions (no tFUS) or how their trial times were made. The same measures are shown for them, all units pooled; this decides no claim.

In [13]:
b = units[units.baseline]
print("Trials per BASELINE recording:", meta[meta.baseline].set_index("identifier").n_trials.to_dict())
print("Units with >2,000 spikes in BASELINE recordings, by recording:", b.groupby("identifier").size().to_dict())
with pd.option_context("display.float_format", "{:.3g}".format):
    print(pd.DataFrame(window_tests(b, "BASELINE")).to_string(index=False))
print("\nBASELINE units with >2,000 spikes, by subject:", b.groupby("subject").size().to_dict(),
      "; total", len(b), "(paper's sham n: 31 S1 RSUs + 53 S1 FSUs = 84)")
Trials per BASELINE recording: {'BH458_BASELINE': 500, 'BH460_BASELINE': 500, 'BH461_BASELINE': 500}
Units with >2,000 spikes in BASELINE recordings, by recording: {'BH458_BASELINE': 69, 'BH460_BASELINE': 16, 'BH461_BASELINE': 34}
     set  window   n  mean_N1  sem_N1  mean_N2  W2_up_p  W1_up_p  W2_down_p  W1_down_p
BASELINE    0-67 119     98.1   0.688     63.7        1        1     0.0351          1
BASELINE  67-134 119      103   0.911     60.8   0.0164        1          1          1
BASELINE 134-201 119      106    1.09     62.4 5.96e-09    0.533          1          1
BASELINE 201-268 119      106    1.36     59.3 1.19e-05    0.943          1          1

BASELINE units with >2,000 spikes, by subject: {'BH458': 69, 'BH460': 16, 'BH461': 34} ; total 119 (paper's sham n: 31 S1 RSUs + 53 S1 FSUs = 84)

Result: BASELINE recordings. The three BASELINE recordings (BH458, BH460, BH461) hold 119 units over the 2,000 spike threshold, against the paper's 84 S1 units in sham (31 RSUs and 53 FSUs). Under W2 they show an increase over the pre-trial window at 134 to 201 ms (corrected p 6.0e-09) and 201 to 268 ms (p 1.2e-05), not under W1 (p 0.53 and 0.94), with mean N1 of 106%. In the time course (trial clock section) their rate at 134 to 335 ms is 104 to 106% of the unit mean against 99% in the pre-window, a size that the BASELINE trace also reaches at times with no trial event (108% at -335 ms, 106% at -201 ms), and there is no rise at 0 to 67 ms. If these were the paper's sham sessions the W2 result at 134 to 268 ms would weigh against fig5-sham-no-response, which tests the 201 to 268 ms window and both tails. Nothing in the files says they are sham sessions (no field describes them, data check), their unit count does not match the paper's sham n, and the claim needs cell types and POm units, so it stays blocked.

Summary¶

Claim Outcome
all 31 claims blocked, data_missing

Every claim is about RSUs and FSUs, which the paper separates by waveform shape; the dataset holds spike times only. Most claims also need POm or auditory cortex units, CG1-targeted or sham recordings, or duty cycles that are not in the dataset (the condition table above). The paper's own data statement names a different dandiset, 000946, whose subjects (20 rats, BH494 to BH594, unit-counts section) are not those of 000637 (BH401 to BH477). No recording of this dataset could be linked to a recording the paper analysed, and the dataset has about 3.6 times as many S1 units over 2,000 spikes in the PRF series at 200 us as the paper's S1 RSUs and FSUs together (4,642 against 1,297; unit-counts section), so these are not, or not only, the units the paper analysed.

What was found about the dataset:

  • The units have no waveforms, electrodes or quality metrics, and all electrode locations are "unknown". The full HDF5 tree of every file was listed: nothing outside spike times, trial times, the electrode table and the subject is stored, and no field names a target, a sham session or a cell type.
  • Trials per recording differ from the 500 the description states (counts in the data check).
  • The description swaps pulse duration and ultrasound duration, and does not list the 600 us recording or the 3.5, 5.5 and 7 V recordings; conditions are coded only in the identifier.
  • Three recordings named BASELINE are not described.
  • Trial starts share a clock with the spikes and behave as tFUS onsets: at 1500 to 4500 Hz, 74% of units fire more in 0 to 67 ms than in the 67 ms before, against 43% with onsets shifted by 1.2 s, and the delayed response is larger at the true onsets in 23 of 25 rats (trial clock section). The elevated rate after the onset lasts about 1 s, so the shifted onsets are a contrast, not a null.
  • Pooling all S1 units (a variant, not a claim), every PRF at 200 us raises the 0 to 67 ms rate over baseline under W2, and 300 to 4500 Hz raise the 201 to 268 ms rate; 30 Hz does not raise it.
  • The three BASELINE recordings (119 units over the threshold, against the paper's 84 S1 units in sham) show an increase over the pre-trial window at 134 to 268 ms under W2 (corrected p 6e-09 and 1.2e-05) and not under W1, of about 5% of the unit mean, a size their trace also reaches at times with no trial event. Nothing says they are sham sessions; if they were, the W2 result would weigh against fig5-sham-no-response.
  • No unit was found deposited twice.

About the claims file: reading N2 of Method S is biased on units with no response (19% at 3 Hz, 116% at 30 Hz, on Poisson units), and the dataset's units have mean N2 of 60 to 86% in all windows but one; under reading R2 it would make claims fail for an artefactual reason. This is asked of a reader in questions/000637.yaml.

How the analysis differs from the paper's: no claim was tested. The pooled-unit variants mix RSUs and FSUs, which the paper does not do, and are shown only to describe the data.