Reproducing Ling et al. (2023/2024) from dandiset 001170¶
Paper. Ling D, Moss EH, Smith CL, Deng F, Kroeger R, Reimer JL, Raman B, Arenkiel BR. Conserved neural dynamics and computations across species in olfaction. bioRxiv doi:10.1101/2023.04.24.538157. This attempt used the full text of version 2 of the preprint (2024-02-23), which is open access (CC BY-ND); its Table S1 and key resources table are images that could not be fetched. The paper was later published in iScience (2026, doi:10.1016/j.isci.2026.115373) with a new title; the claims tested here are those of the preprint, which the dandiset names.
Dandiset. DANDI:001170, version 0.260528.1016. The paper cites no DANDI version (its data statement promises Figshare). This is the most recent published version; the only other one, 0.251111.1532, holds the same 12 assets (same asset ids, sizes and digests, shown in the data check) and differs only in metadata (anatomy terms, keywords, funders, affiliations and a longer description), none of which the analysis uses.
What the dandiset holds. Two-photon imaging of the dorsal olfactory bulb of awake mice: raw fluorescence of glomerular ROIs, the raw movies, a trials table of odor valve openings, and treadmill (and in one file respiration) traces. Seven sessions from four mice follow the paper's block protocol (Figures 2, 4, 6, 7); five sessions from five mice (one of them, sub-134, also in the block sessions) follow a random-order protocol with 11 odors. The locust recordings of the paper are not in the dandiset; they are on Figshare with the paper's code.
Outside inputs. The paper's Figshare record
(doi:10.6084/m9.figshare.31435543.v2): its
MATLAB scripts (read for method choices), the deposited mouse array data1_mouse_glom198.mat
(used only by the link claim), and the locust files data2_locust_PN66.mat,
data4_locust_PN85.mat and data5_locust_PN85.mat (used by the locust claims, data elsewhere).
Each file is checked against the md5 the record gives.
Random numbers (random vectors and label permutations) come from a generator seeded per claim.
| Claim | Data | Status |
|---|---|---|
| mouse-sessions | dandiset | tested |
| figshare-trace-link | with inputs | tested |
| fig2h-onoff-angle-mouse | dandiset | tested |
| fig2f-onoff-anticorrelated-mouse | dandiset | tested |
| fig4b-reference-odor-mouse | dandiset | tested |
| fig4d-clustering-mouse | dandiset | tested |
| fig6e-contrast-msc | dandiset | tested |
| fig6e-contrast-ace | dandiset | blocked: no sequence ends in acetophenone |
| fig7d-persistent-anticorrelation | dandiset | tested |
| fig7h-iti-knn-mouse | dandiset | tested |
| fig7j-iti-timebin-mouse | dandiset | tested |
| figs7e-ibi-chance | dandiset | tested |
| figs8e-random-order-decline | dandiset | tested |
| figs2-pid-decay | dandiset | blocked: no PID recordings |
| fig1g-onoff-angle-locust | elsewhere (Figshare) | tested |
| figs4a-reference-odor-locust | elsewhere (Figshare) | tested |
| fig5e-contrast-locust | elsewhere (Figshare) | tested |
| fig7f-iti-knn-locust | elsewhere (Figshare) | tested |
| fig7i-iti-timebin-locust | elsewhere (Figshare) | tested |
import os
os.environ.setdefault("OMP_NUM_THREADS", "4")
import time
import warnings
from concurrent.futures import ThreadPoolExecutor
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
from scipy import stats
from scipy.cluster.hierarchy import fcluster, linkage
from scipy.spatial.distance import pdist
import olfdyn
pd.set_option("display.width", 220)
pd.set_option("display.max_columns", 40)
pd.set_option("display.max_rows", 200)
warnings.filterwarnings("ignore", category=RuntimeWarning)
started = time.time()
urls = olfdyn.asset_urls()
with ThreadPoolExecutor(3) as pool:
sessions = list(pool.map(lambda p: olfdyn.reduce_session(p, urls[p]), sorted(urls)))
for s in sessions:
s["name"] = f"sub-{s['subject']} ses-{s['session']}"
block = [s for s in sessions if s["protocol"] == "block"]
random_sessions = [s for s in sessions if s["protocol"] == "random"]
print(f"loaded {len(sessions)} files from dandiset {olfdyn.DANDISET_ID} version {olfdyn.DANDISET_VERSION} "
f"in {time.time() - started:.0f} s")
loaded 12 files from dandiset 001170 version 0.260528.1016 in 72 s
Data check¶
Each check of the metadata against the values is stated under its output.
from dandi.dandiapi import DandiAPIClient
with DandiAPIClient() as client:
versions = [v.identifier for v in client.get_dandiset(olfdyn.DANDISET_ID).get_versions()]
print("Versions on the archive:", versions)
tables = {v: olfdyn.asset_table(v) for v in versions if v != "draft"}
first, last = tables["0.251111.1532"], tables[olfdyn.DANDISET_VERSION]
print("Assets per version:", {v: len(t) for v, t in tables.items()})
print("0.251111.1532 and", olfdyn.DANDISET_VERSION, "hold the same paths, asset ids, sizes and digests:",
first.equals(last))
print(f"Total size: {last['size'].sum() / 1e9:.1f} GB")
Versions on the archive: ['draft', '0.251111.1532', '0.260528.1016']
Assets per version: {'0.251111.1532': 12, '0.260528.1016': 12}
0.251111.1532 and 0.260528.1016 hold the same paths, asset ids, sizes and digests: True
Total size: 352.1 GB
info = pd.DataFrame([{
"session": s["name"], "sex": s["sex"], "dob": s["dob"], "date": s["start"][:10], "protocol": s["protocol"],
"rois": s["n_rois"], "rois_per_fov": s["rois_per_fov"], "frames": s["n_frames"],
"frame_hz": round(s["frame_hz"], 2), "unit": s["unit"], "indicator": s["indicator"],
"location": s["location"], "same_t_3fov": s["same_timestamps"], "behavior": ",".join(s["behavior"]),
"description": s["subject_description"]} for s in sessions])
print(info.to_string(index=False))
print("\nSessions per protocol:", info.protocol.value_counts().to_dict())
print("Mice per protocol:", info.assign(m=info.session.str.split().str[0]).groupby("protocol").m.nunique().to_dict())
print("Frame rates (Hz): min", info.frame_hz.min(), "max", info.frame_hz.max(), "(Methods: 15-18 Hz; code: 16)")
print("Session start times:", sorted({s['start'][11:] for s in sessions}))
session sex dob date protocol rois rois_per_fov frames frame_hz unit indicator location same_t_3fov behavior description
sub-124 ses-3 M 2020-11-28 2021-03-22 random 559 [182, 237, 140] 73610 15.28 n.a. GCaMP6f OB True respiration,treadmill_position ArenkielMouse
sub-125 ses-1 F 2020-11-29 2021-03-16 random 325 [105, 142, 78] 147000 15.10 n.a. GCaMP6f OB True treadmill_position ArenkielMouse
sub-126 ses-1 M 2021-02-10 2021-04-09 random 323 [96, 132, 95] 133000 13.73 n.a. GCaMP6f OB True treadmill_position ArenkielMouse
sub-134 ses-22 M 2021-07-15 2022-07-21 block 73 [17, 48, 8] 182000 17.83 n.a. GCaMP6f OB True treadmill_position ArenkielMouse
sub-134 ses-6 M 2021-07-15 2022-01-11 random 306 [99, 148, 59] 75000 15.29 n.a. GCaMP6f OB True treadmill_position ArenkielMouse
sub-141 ses-1 M 2022-02-24 2022-05-12 block 210 [56, 102, 52] 150000 14.95 n.a. GCaMP6f OB True treadmill_position ArenkielMouse
sub-141 ses-5 M 2022-02-24 2022-06-24 block 153 [71, 50, 32] 147000 14.40 n.a. GCaMP6f OB True treadmill_position ArenkielMouse
sub-141 ses-8 M 2022-02-24 2022-07-08 block 198 [80, 87, 31] 167000 15.50 n.a. GCaMP6f OB True treadmill_position ArenkielMouse
sub-148 ses-2 F 2022-07-01 2023-02-28 block 286 [95, 106, 85] 165000 16.38 n.a. GCaMP6f OB True treadmill_position Arenkiel Mouse
sub-148 ses-3 F 2022-07-01 2023-03-01 block 318 [106, 127, 85] 160000 15.72 n.a. GCaMP6f OB True treadmill_position Arenkiel Mouse
sub-149 ses-1 F 2022-07-01 2023-02-08 block 219 [63, 100, 56] 155000 15.29 n.a. GCaMP6f OB True treadmill_position Arenkiel Mouse
sub-1716 ses-1 M 2020-03-11 2020-07-31 random 403 [143, 91, 169] 140000 14.24 n.a. GCaMP6f OB True treadmill_position arenkiel mouse
Sessions per protocol: {'block': 7, 'random': 5}
Mice per protocol: {'block': 4, 'random': 5}
Frame rates (Hz): min 13.73 max 17.83 (Methods: 15-18 Hz; code: 16)
Session start times: ['00:00:00-05:00', '00:00:00-06:00']
Twelve sessions from eight mice. The fluorescence is raw F (unit "n.a."), not the dF/F the Methods describe, and no description says how it was extracted. The indicator (GCaMP6f) and the 920 nm excitation agree with the Methods; the genotype (Thy1-GCaMP6f) and the ages (two to four months) are not stated, though dates of birth are. Every session starts at 00:00, so the time of day is missing. The three fields of view of a session share one timestamp array. Frame rates are 13.7 to 17.8 Hz, against the Methods' 15 to 18 Hz and the code's nominal 16.
rows = []
for s in sessions:
tp = s["two_photon"]
fl_dt = 1 / s["frame_hz"]
for k, (n, dt) in tp.items():
rows.append({"session": s["name"], "series": k.replace("two_photon_series_", ""), "movie_frames": n,
"trace_frames": s["n_frames"], "ratio": round(s["n_frames"] / n, 3),
"movie_dt_ms": round(1000 * dt, 1), "trace_dt_ms": round(1000 * fl_dt, 1)})
tp = pd.DataFrame(rows)
print(tp[tp.series.str.startswith("FOV1")].to_string(index=False))
print("\nSessions whose raw movie holds fewer frames than the traces:",
tp[tp.ratio != 1].session.unique().tolist())
session series movie_frames trace_frames ratio movie_dt_ms trace_dt_ms sub-124 ses-3 FOV1_channel1 73610 73610 1.0 65.6 65.4 sub-125 ses-1 FOV1_channel1 147000 147000 1.0 66.4 66.2 sub-126 ses-1 FOV1_channel1 133000 133000 1.0 73.0 72.8 sub-134 ses-22 FOV1_channel1 91000 182000 2.0 112.4 56.1 sub-134 ses-6 FOV1_channel1 75000 75000 1.0 65.5 65.4 sub-141 ses-1 FOV1_channel1 150000 150000 1.0 67.0 66.9 sub-141 ses-5 FOV1_channel1 73500 147000 2.0 139.1 69.4 sub-141 ses-8 FOV1_channel1 83500 167000 2.0 129.2 64.5 sub-148 ses-2 FOV1_channel1 165000 165000 1.0 61.2 61.0 sub-148 ses-3 FOV1_channel1 160000 160000 1.0 63.8 63.6 sub-149 ses-1 FOV1_channel1 155000 155000 1.0 65.4 65.4 sub-1716 ses-1 FOV1_channel1 140000 140000 1.0 70.4 70.2 Sessions whose raw movie holds fewer frames than the traces: ['sub-134 ses-22', 'sub-141 ses-5', 'sub-141 ses-8']
In three sessions (sub-134 ses-22, sub-141 ses-5, sub-141 ses-8) the raw TwoPhotonSeries holds half the frames of the fluorescence traces over the same span, at twice the frame interval. The conversion notes in the conversion repository give 182000 frames at 17.8 Hz for sub-134 ses-22, the count of its traces, so the movies of these three files are missing every other frame (or were written at half rate). The analysis uses the traces only.
rows = []
for s in sessions:
tr, p = s["trials"], s["pres"]
dur = tr.stop_time - tr.start_time
rows.append({"session": s["name"], "protocol": s["protocol"], "rows": len(tr),
"empty_rows": int((tr.odorant == olfdyn.EMPTY).sum()),
"row_s": f"{dur.min():.2f}-{dur.max():.2f}", "presentations": len(p),
"pulse_rows": sorted(p.rows.unique().tolist()),
"odors": len(p.odor.unique()), "concentrations": sorted(tr.concentration.unique().tolist()),
"conc_of_empty": sorted(tr.concentration[tr.odorant == olfdyn.EMPTY].unique().tolist()),
"first_trial_s": round(tr.start_time.min() - s["clock"]["t_first"], 1),
"after_last_s": round(s["clock"]["t_last"] - tr.stop_time.max(), 1)})
print(pd.DataFrame(rows).to_string(index=False))
names = sorted({n for s in sessions for n in s["trials"].odorant.unique()})
print("\nOdorant names in the files:", names)
session protocol rows empty_rows row_s presentations pulse_rows odors concentrations conc_of_empty first_trial_s after_last_s sub-124 ses-3 random 241 15 1.00-1.01 226 [1] 11 [0.0, 0.01, 0.05, 0.1, 0.6, 0.7, 0.8, 0.9] [0.0] 34.1 -8.9 sub-125 ses-1 random 480 53 1.00-1.01 427 [1] 11 [0.0, 0.01, 0.05, 0.1, 0.6, 0.7, 0.8, 0.9] [0.0] 39.6 136.5 sub-126 ses-1 random 480 53 1.00-1.01 427 [1] 11 [0.0, 0.01, 0.05, 0.1, 0.6, 0.7, 0.8, 0.9] [0.0] 30.4 90.9 sub-134 ses-22 block 500 0 1.00-1.01 100 [5] 4 [0.1] [] 19.1 591.5 sub-134 ses-6 random 240 15 1.00-1.01 225 [1] 11 [0.0, 0.01] [0.0] 15.2 112.0 sub-141 ses-1 block 500 0 1.00-1.01 100 [5] 4 [0.1] [] 23.8 413.9 sub-141 ses-5 block 9900 9400 1.00-1.01 100 [5] 4 [0.1] [0.1] 14.4 237.7 sub-141 ses-8 block 9900 9400 1.00-1.10 100 [5] 4 [0.0, 0.1] [0.0] 19.7 803.8 sub-148 ses-2 block 9900 9400 1.00-1.10 100 [5] 4 [0.0, 0.1] [0.0] 24.6 101.8 sub-148 ses-3 block 9900 9400 1.00-1.10 100 [5] 4 [0.0, 0.1] [0.0] 34.2 199.2 sub-149 ses-1 block 9900 9400 1.00-1.10 100 [5] 4 [0.0, 0.1] [0.0] 13.3 161.6 sub-1716 ses-1 random 480 53 1.00-1.01 427 [1] 11 [0.0, 0.01] [0.0] 26.8 246.3 Odorant names in the files: ['1,4-Cineole', '1,4-cineole', '2-Ethylphenol', '2-ethylphenol', '2-methypyrazine', 'Acetophenone', 'Allyl Sulfide', 'Anisole', 'Empty', 'Eugenol', 'Isoeugenol', 'Methyl Salicylate', 'anisaldehyde', 'anisole', 'trans-Cinnamaldehyde', 'trans-cinnamaldehyde']
Each trials-table row is one second of valve opening: a 5-s presentation is five rows, and in some sessions the inter-trial seconds are rows named "Empty" (sub-141 ses-5 gives these rows a concentration of 0.1, the other sessions 0). Two block sessions (sub-134 ses-22, sub-141 ses-1) have no Empty rows. Presentations are rebuilt here by joining consecutive rows of one odorant. The concentration has no unit. Odorant names are spelled differently across sessions (for example "Anisole" and "anisole"); "2-methypyrazine" is a misspelling. Names of the random-order odors differ from the paper's list: the files have 2-ethylphenol and anisaldehyde where the Methods list 2-ethylpyrizine (2ep) and isoamyl acetate (asa).
for s in block:
p = s["pres"]
summary = []
for b, g in p.groupby("block"):
kinds = g.kind.unique().tolist()
odors = "/".join(g[g.kind != "seq_second"].odor.unique())
second = "/".join(g[g.kind == "seq_second"].odor.unique())
pulse = (g.offset - g.onset).median()
isi = np.median(np.diff(g[g.kind != "seq_second"].onset)) if len(g) > 1 else np.nan
pair = (g[g.kind == "seq_second"].onset.to_numpy() - g[g.kind == "seq_first"].onset.to_numpy()).mean() \
if "seq_first" in kinds else np.nan
summary.append(f"{'sol' if kinds == ['solitary'] else 'seq'}:{odors}{('->' + second) if second else ''}"
f" n={int((g.kind != 'seq_second').sum())} pulse={pulse:.2f}s isi={isi:.1f}s"
+ (f" pair_lag={pair:.2f}s" if pair == pair else ""))
gaps = np.diff([g.onset.min() for _, g in p.groupby("block")])
ends = [g.onset.max() for _, g in p.groupby("block")]
ibis = [p[p.block == b + 1].onset.min() - ends[b] for b in range(len(ends) - 1)]
print(s["name"], "|", "; ".join(summary))
print(" time from a block's last onset to the next block's first onset (s):", np.round(ibis, 0).tolist())
sub-134 ses-22 | sol:msc n=10 pulse=5.03s isi=60.4s; sol:ace n=10 pulse=5.03s isi=60.4s; sol:als n=10 pulse=5.03s isi=60.4s; sol:eug n=10 pulse=5.03s isi=60.4s; seq:ace->msc n=10 pulse=5.03s isi=60.4s pair_lag=6.05s; seq:als->msc n=10 pulse=5.03s isi=60.4s pair_lag=6.04s; seq:eug->msc n=10 pulse=5.03s isi=60.4s pair_lag=6.04s
time from a block's last onset to the next block's first onset (s): [967.0, 967.0, 967.0, 967.0, 961.0, 961.0]
sub-141 ses-1 | sol:msc n=10 pulse=5.03s isi=60.4s; sol:ace n=10 pulse=5.03s isi=60.4s; sol:als n=10 pulse=5.03s isi=60.4s; sol:eug n=10 pulse=5.03s isi=60.4s; seq:ace->msc n=10 pulse=5.03s isi=60.4s pair_lag=6.04s; seq:als->msc n=10 pulse=5.03s isi=60.4s pair_lag=6.04s; seq:eug->msc n=10 pulse=5.03s isi=60.4s pair_lag=6.04s
time from a block's last onset to the next block's first onset (s): [967.0, 967.0, 967.0, 967.0, 961.0, 961.0]
sub-141 ses-5 | sol:msc n=10 pulse=5.03s isi=60.4s; sol:ace n=10 pulse=5.03s isi=60.4s; sol:als n=10 pulse=5.03s isi=60.4s; sol:eug n=10 pulse=5.03s isi=60.4s; seq:ace->msc n=10 pulse=5.03s isi=60.4s pair_lag=6.05s; seq:als->msc n=10 pulse=5.03s isi=60.4s pair_lag=6.04s; seq:eug->msc n=10 pulse=5.03s isi=60.4s pair_lag=6.04s
time from a block's last onset to the next block's first onset (s): [967.0, 967.0, 967.0, 967.0, 961.0, 961.0]
sub-141 ses-8 | sol:msc n=10 pulse=5.03s isi=60.4s; sol:ace n=10 pulse=5.03s isi=60.4s; sol:als n=10 pulse=5.03s isi=60.4s; sol:eug n=10 pulse=5.03s isi=60.4s; seq:ace->msc n=10 pulse=5.03s isi=60.4s pair_lag=6.05s; seq:als->msc n=10 pulse=5.03s isi=60.4s pair_lag=6.04s; seq:eug->msc n=10 pulse=5.03s isi=60.4s pair_lag=6.04s
time from a block's last onset to the next block's first onset (s): [967.0, 967.0, 967.0, 967.0, 961.0, 961.0]
sub-148 ses-2 | seq:ace->msc n=10 pulse=5.03s isi=60.4s pair_lag=6.04s; seq:als->msc n=10 pulse=5.03s isi=60.5s pair_lag=6.04s; seq:eug->msc n=10 pulse=5.03s isi=60.4s pair_lag=6.05s; sol:msc n=10 pulse=5.04s isi=60.4s; sol:ace n=10 pulse=5.03s isi=60.4s; sol:als n=10 pulse=5.03s isi=60.4s; sol:eug n=10 pulse=5.03s isi=60.5s
time from a block's last onset to the next block's first onset (s): [961.0, 961.0, 913.0, 967.0, 967.0, 967.0]
sub-148 ses-3 | seq:ace->msc n=10 pulse=5.03s isi=60.4s pair_lag=6.04s; seq:als->msc n=10 pulse=5.03s isi=60.4s pair_lag=6.04s; seq:eug->msc n=10 pulse=5.03s isi=60.4s pair_lag=6.05s; sol:msc n=10 pulse=5.03s isi=60.4s; sol:ace n=10 pulse=5.03s isi=60.4s; sol:als n=10 pulse=5.04s isi=60.4s; sol:eug n=10 pulse=5.04s isi=60.5s
time from a block's last onset to the next block's first onset (s): [961.0, 961.0, 913.0, 967.0, 967.0, 967.0]
sub-149 ses-1 | seq:ace->msc n=10 pulse=5.03s isi=60.4s pair_lag=6.04s; seq:als->msc n=10 pulse=5.03s isi=60.4s pair_lag=6.04s; seq:eug->msc n=10 pulse=5.03s isi=60.4s pair_lag=6.04s; sol:msc n=10 pulse=5.03s isi=60.4s; sol:ace n=10 pulse=5.04s isi=60.4s; sol:als n=10 pulse=5.03s isi=60.4s; sol:eug n=10 pulse=5.03s isi=60.4s
time from a block's last onset to the next block's first onset (s): [961.0, 961.0, 913.0, 967.0, 967.0, 967.0]
All seven block sessions hold four solitary blocks (methyl salicylate, acetophenone, allyl sulfide, eugenol; ten 5-s pulses every 60 s, so 55 s between pulses) and three sequential blocks in which acetophenone, allyl sulfide or eugenol is followed 6 s after its onset (1 s after its offset) by methyl salicylate, ten pairs every 60 s. Blocks are about 15 min apart (the last onset of a block to the next block's first onset is 961 to 967 s, and 913 s between the last sequential and the first solitary block in sub-148 and sub-149). This agrees with the Methods. No block presents acetophenone, or any odor other than methyl salicylate, as the second odor of a pair. The order of blocks differs between mice: sub-134 and sub-141 start with the solitary blocks, sub-148 and sub-149 with the sequential ones.
for s in random_sessions:
p = s["pres"]
isi = np.diff(p.onset)
print(f"{s['name']}: {len(p)} presentations of {p.odor.nunique()} odorants incl. Empty, pulse "
f"{(p.offset - p.onset).median():.2f} s, onset interval median {np.median(isi):.2f} s "
f"(range {isi.min():.2f}-{isi.max():.2f}); presentations per odor (Empty excluded) "
f"{p[p.odorant != olfdyn.EMPTY].odor.value_counts().min()}-{p[p.odorant != olfdyn.EMPTY].odor.value_counts().max()}")
print(" concentrations by odor:",
{o: sorted(g.concentration.unique().tolist()) for o, g in p[p.odorant != olfdyn.EMPTY].groupby("odor")})
seqs = {s["name"]: tuple(s["pres"].odor.str.lower()) for s in random_sessions}
base = seqs[random_sessions[0]["name"]]
print("\nRandom-order sessions whose first 200 presentations match", random_sessions[0]["name"], ":",
[n for n, q in seqs.items() if q[:200] == base[:200]])
sub-124 ses-3: 226 presentations of 11 odorants incl. Empty, pulse 1.00 s, onset interval median 20.01 s (range 20.00-60.02); presentations per odor (Empty excluded) 13-26
concentrations by odor: {'2ep': [0.0], 'ace': [0.05], 'als': [0.8], 'ani': [0.7], 'asa': [0.0], 'cin': [0.01], 'eug': [0.9], 'ieg': [0.9], 'mpz': [0.01], 'msc': [0.6], 'tcm': [0.1]}
sub-125 ses-1: 427 presentations of 11 odorants incl. Empty, pulse 1.00 s, onset interval median 20.01 s (range 20.00-100.03); presentations per odor (Empty excluded) 13-67
concentrations by odor: {'2ep': [0.0], 'ace': [0.05], 'als': [0.8], 'ani': [0.7], 'asa': [0.0], 'cin': [0.01], 'eug': [0.9], 'ieg': [0.9], 'mpz': [0.01], 'msc': [0.6], 'tcm': [0.1]}
sub-126 ses-1: 427 presentations of 11 odorants incl. Empty, pulse 1.00 s, onset interval median 20.01 s (range 20.00-100.03); presentations per odor (Empty excluded) 13-67
concentrations by odor: {'2ep': [0.0], 'ace': [0.05], 'als': [0.8], 'ani': [0.7], 'asa': [0.0], 'cin': [0.01], 'eug': [0.9], 'ieg': [0.9], 'mpz': [0.01], 'msc': [0.6], 'tcm': [0.1]}
sub-134 ses-6: 225 presentations of 11 odorants incl. Empty, pulse 1.00 s, onset interval median 20.01 s (range 20.00-60.02); presentations per odor (Empty excluded) 13-26
concentrations by odor: {'2ep': [0.01], 'ace': [0.01], 'als': [0.01], 'ani': [0.01], 'asa': [0.01], 'cin': [0.01], 'eug': [0.01], 'ieg': [0.01], 'mpz': [0.01], 'msc': [0.01], 'tcm': [0.01]}
sub-1716 ses-1: 427 presentations of 11 odorants incl. Empty, pulse 1.00 s, onset interval median 20.01 s (range 20.00-100.03); presentations per odor (Empty excluded) 13-67
concentrations by odor: {'2ep': [0.01], 'ace': [0.01], 'als': [0.01], 'ani': [0.01], 'asa': [0.01], 'cin': [0.01], 'eug': [0.01], 'ieg': [0.01], 'mpz': [0.01], 'msc': [0.01], 'tcm': [0.01]}
Random-order sessions whose first 200 presentations match sub-124 ses-3 : ['sub-124 ses-3', 'sub-125 ses-1', 'sub-126 ses-1', 'sub-134 ses-6', 'sub-1716 ses-1']
The five random-order sessions present 1-s pulses every 20 s, as the Methods' random-order experiment does, and share one fixed odor sequence. In three of them (sub-124, sub-125, sub-126) each odor has its own concentration value, from 0.01 to 0.9, and two odors (2-ethylphenol and anisaldehyde) have 0.0, the value given to the blank "Empty"; the Methods give one dilution for all odors. The other two sessions give 0.01 to every odor. None of the five has the 286 glomeruli the Methods give for Figure S8D. In sub-124 ses-3 the last trial ends 8.9 s after the last imaging frame.
# Recordings deposited twice: frame counts, ROI counts and the first frames of the traces.
keys = pd.DataFrame([{"session": s["name"], "frames": s["n_frames"], "rois": s["n_rois"],
"first_frame_mean": float(np.round(s["F_median"].mean(), 3))} for s in sessions])
dup_frames = keys[keys.duplicated("frames", keep=False)]
print("Sessions sharing a frame count:\n", dup_frames.to_string(index=False))
pairs = []
for i, a in enumerate(sessions):
for b in sessions[i + 1:]:
if a["n_rois"] == b["n_rois"] or np.isclose(a["F_median"].mean(), b["F_median"].mean()):
pairs.append((a["name"], b["name"]))
print("Pairs sharing a ROI count or the same mean of median F:", pairs)
same_trials = [(a["name"], b["name"]) for i, a in enumerate(sessions) for b in sessions[i + 1:]
if len(a["trials"]) == len(b["trials"]) and np.allclose(a["trials"].start_time, b["trials"].start_time)]
print("Pairs with identical trial start times:", same_trials)
Sessions sharing a frame count:
session frames rois first_frame_mean
sub-125 ses-1 147000 325 87.917999
sub-141 ses-5 147000 153 544.379028
Pairs sharing a ROI count or the same mean of median F: []
Pairs with identical trial start times: []
Two sessions share a frame count (sub-125 ses-1 and sub-141 ses-5, 147000), but they differ in ROI count and protocol. No two sessions share ROI counts or trial times, so no recording appears to be deposited twice. Sessions of one mouse (sub-134, sub-141, sub-148) are analysed separately; nothing says that their ROIs are the same glomeruli.
# Clock: the traces' timestamps and the trials table are both on the odor clock (conversion code).
# Block sessions: latency of the population response after each onset, from the trial-averaged
# absolute change of F relative to the 2 s before onset.
rows = []
for s in block:
x = olfdyn.dff_pre(s["sol"])
curve = np.abs(x.mean(axis=0)).mean(axis=1) # frames
pre = curve[:32]
thr = pre.mean() + 5 * pre.std()
after = np.flatnonzero(curve[32:] > thr)
rows.append({"session": s["name"], "latency_s": (after[0] / 16) if len(after) else np.nan,
"peak_s": float(np.argmax(curve[32:112]) / 16), "pre_level": float(pre.mean()),
"peak_level": float(curve[32:112].max())})
print(pd.DataFrame(rows).round(3).to_string(index=False))
# Random-order sessions: decode the odor from 0.25 to 2 s after onset, with the true labels and with
# the labels of the presentation before (a contrast between conditions).
rows = []
for s in random_sessions:
q = s["rnd_pres"]
x = olfdyn.dff_pre(s["rnd"])[:, 36:64].mean(axis=1)
nb = olfdyn.knn_neighbours(x, np.arange(len(x)), 5, "cosine")
true = q.odor.to_numpy()
prev = q.previous.fillna("none").to_numpy()
rows.append({"session": s["name"], "n": len(q), "acc_true": (olfdyn.knn_predict(nb, true) == true).mean(),
"acc_previous_label": (olfdyn.knn_predict(nb, prev) == prev).mean(), "chance": 1 / 11})
print(pd.DataFrame(rows).round(3).to_string(index=False))
session latency_s peak_s pre_level peak_level sub-134 ses-22 0.312 4.188 0.013 0.079 sub-141 ses-1 0.250 0.688 0.011 0.233 sub-141 ses-5 0.250 4.750 0.011 0.081 sub-141 ses-8 0.312 4.625 0.013 0.159 sub-148 ses-2 0.250 0.562 0.007 0.046 sub-148 ses-3 0.188 0.500 0.007 0.076 sub-149 ses-1 0.250 0.875 0.011 0.068
session n acc_true acc_previous_label chance sub-124 ses-3 110 0.255 0.064 0.091 sub-125 ses-1 110 0.764 0.064 0.091 sub-126 ses-1 110 0.545 0.055 0.091 sub-134 ses-6 110 0.655 0.073 0.091 sub-1716 ses-1 110 0.836 0.064 0.091
In every block session the population response rises 0.19 to 0.31 s after the valve onset and peaks within the 5 s pulse (0.5 to 4.75 s). In the random-order sessions the odor is decoded from the 2 s after onset at 25.5% to 83.6% with the true labels and at 5.5% to 7.3% with the label of the previous presentation (chance 9.1%). This is consistent with the traces and the trials table sharing a clock; a constant offset shorter than about half a second cannot be excluded from these outputs.
# Grid points of the analysis windows that fall before the first or after the last imaging frame
# (left as NaN by olfdyn.windows, not extrapolated).
rows = []
for s in sessions:
row = {"session": s["name"], "protocol": s["protocol"]}
for key in ["sol", "seq", "rnd"]:
if key in s:
bad = np.isnan(s[key][..., 0])
row[key + "_trials_cut"] = int(bad.any(axis=1).sum())
row[key + "_points_cut"] = int(bad.sum())
if "pseudo" in s:
row["pseudo_points_cut"] = int(sum(np.isnan(v[..., 0]).sum() for v in s["pseudo"].values()))
rows.append(row)
print(pd.DataFrame(rows).fillna("-").to_string(index=False))
session protocol rnd_trials_cut rnd_points_cut sol_trials_cut sol_points_cut seq_trials_cut seq_points_cut pseudo_points_cut sub-124 ses-3 random 0.0 0.0 - - - - - sub-125 ses-1 random 0.0 0.0 - - - - - sub-126 ses-1 random 0.0 0.0 - - - - - sub-134 ses-22 block - - 0.0 0.0 0.0 0.0 0.0 sub-134 ses-6 random 0.0 0.0 - - - - - sub-141 ses-1 block - - 0.0 0.0 0.0 0.0 0.0 sub-141 ses-5 block - - 0.0 0.0 0.0 0.0 0.0 sub-141 ses-8 block - - 0.0 0.0 0.0 0.0 0.0 sub-148 ses-2 block - - 0.0 0.0 0.0 0.0 0.0 sub-148 ses-3 block - - 0.0 0.0 0.0 0.0 0.0 sub-149 ses-1 block - - 0.0 0.0 0.0 0.0 0.0 sub-1716 ses-1 random 0.0 0.0 - - - - -
No grid point of any analysis window falls before the first or after the last imaging frame, so no value is extrapolated. The trials of sub-124 ses-3 that run past the end of imaging are not among the first ten presentations of each odor used here.
mouse-sessions¶
The mouse block-protocol data (four odors, 5 s pulses, 55 s between trials, blocks of ten, 15 min between blocks) come from seven imaging sessions in four mice, with a median of 215 glomeruli per imaging session.
n_block = len(block)
mice = sorted({s["subject"] for s in block})
rois_block = [s["n_rois"] for s in block]
rois_all = [s["n_rois"] for s in sessions]
print("Block-protocol sessions:", n_block, [s["name"] for s in block])
print("Mice:", len(mice), mice)
print("ROIs per block session:", dict(zip([s["name"] for s in block], rois_block)))
print("(a) median over the block sessions:", np.median(rois_block))
print("(b) median over all 12 sessions:", np.median(rois_all), "range", min(rois_all), "-", max(rois_all))
print("Expected 7 sessions, 4 mice, median 215 +/- 21.5")
Block-protocol sessions: 7 ['sub-134 ses-22', 'sub-141 ses-1', 'sub-141 ses-5', 'sub-141 ses-8', 'sub-148 ses-2', 'sub-148 ses-3', 'sub-149 ses-1']
Mice: 4 ['134', '141', '148', '149']
ROIs per block session: {'sub-134 ses-22': 73, 'sub-141 ses-1': 210, 'sub-141 ses-5': 153, 'sub-141 ses-8': 198, 'sub-148 ses-2': 286, 'sub-148 ses-3': 318, 'sub-149 ses-1': 219}
(a) median over the block sessions: 210.0
(b) median over all 12 sessions: 296.0 range 73 - 559
Expected 7 sessions, 4 mice, median 215 +/- 21.5
Result: reproduced. The trials tables give seven block-protocol sessions from four mice (134, 141, 148, 149). Over these sessions the median ROI count is 210, within 215 +/- 21.5 (reading a); over all twelve sessions it is 296 (reading b). The ROI counts of sub-141 ses-8 (198), sub-149 ses-1 (219) and sub-148 ses-3 (318) equal the glomerulus counts the Methods give for Figures 2F and 4C, 6D, and 7, so these sessions are taken to be those figures' sessions. This is an assumption by count: no identifier links the paper's datasets to the dandiset's sessions.
figshare-trace-link¶
The dandiset's raw fluorescence of the 198-glomerulus session is the data the paper deposited for its Figure 2D to 2F and 4C (Figshare data1_mouse_glom198.mat), up to the normalisation, with glomeruli in the same order.
Result: reproduced. With windows counted in frames as the code does, the median within-trial r between the dandiset's raw F of sub-141 ses-8 and the deposited array is 0.968 (5th percentile 0.834; 0.947 to 0.994 by block), and every deposited glomerulus correlates best with the dandiset ROI of the same index. With windows placed in seconds the median r is 0.307: the code counts 16 frames per second in a session imaged at 15.5 Hz, so its windows drift by about 1.8 s over a 912-frame trial. The deposited values are close to F z-scored within each trial window (r = 0.927).
Responses of the block sessions¶
Trial-averaged responses to the four solitary odors in each block session, used by the next claims. Frames are on a 16 Hz grid from 2 s before onset; ON is frames 32 to 112 (0 to 5 s after onset), OFF is frames 112 to 192 (0 to 5 s after offset).
PRE, ON, OFF = slice(0, 32), slice(32, 112), slice(112, 192)
resp = {} # (session, odor) -> trial-averaged dF/F, frames x ROIs (F0 = 2 s before onset)
resp_z = {} # the same from per-trial z-scored F (variant)
resp_zpre = {} # the same z-scored with the mean and SD of the 2 s before onset only (variant)
resp_fr = {} # dF/F on frames counted as the code does, 16 per nominal second (variant)
for s in block:
x = olfdyn.dff_pre(s["sol"])
z = olfdyn.zscore_trials(s["sol"][:, :912])
zp = olfdyn.zscore_pre(s["sol"])
fr = s["sol_frames"].astype(np.float64)
fr = (fr - fr[:, :32].mean(axis=1, keepdims=True)) / fr[:, :32].mean(axis=1, keepdims=True)
for odor in olfdyn.BLOCK_ODORS:
sel = (s["sol_pres"].odor == odor).to_numpy()
resp[s["name"], odor] = x[sel].mean(axis=0).astype(np.float64)
resp_z[s["name"], odor] = z[sel].mean(axis=0).astype(np.float64)
resp_zpre[s["name"], odor] = zp[sel].mean(axis=0).astype(np.float64)
resp_fr[s["name"], odor] = fr[sel].mean(axis=0)
print("NaN values in the trial-averaged responses:", int(sum(np.isnan(v).sum() for v in resp.values())))
print("Trials per session and odor:", sorted({int((s['sol_pres'].odor == o).sum()) for s in block for o in olfdyn.BLOCK_ODORS}))
print("dF/F range over all trial-averaged frames:", round(min(v.min() for v in resp.values()), 3), "to",
round(max(v.max() for v in resp.values()), 3))
print("F0 (2 s pre-onset mean) range over trials and ROIs:",
round(min(s["sol"][:, PRE].mean(axis=1).min() for s in block), 1), "to",
round(max(s["sol"][:, PRE].mean(axis=1).max() for s in block), 1))
NaN values in the trial-averaged responses: 0 Trials per session and odor: [10] dF/F range over all trial-averaged frames: -0.43 to 2.523 F0 (2 s pre-onset mean) range over trials and ROIs: 60.9 to 2295.0
fig2h-onoff-angle-mouse¶
In the mouse OB the angle between glomerular ON and OFF response vectors is greater than 90 degrees for each odor (mean 121.8 +/- 2.98 SEM, 4 odors, 7 sessions), so ON and OFF responses are inverted relative to each other.
rng = np.random.default_rng(2)
rows = []
for s in block:
for odor in olfdyn.BLOCK_ODORS:
r, rz = resp[s["name"], odor], resp_z[s["name"], odor]
med = olfdyn.median_subtracted(r)
main = olfdyn.onoff_angles(med[ON], med[OFF], center=True)
noise = rng.standard_normal((160, r.shape[1]))
rows.append({"session": s["name"], "mouse": s["subject"], "odor": odor, "on_off": main["on_off"],
"on_on": main["on_on"], "off_off": main["off_off"],
"uncentred": olfdyn.onoff_angles(med[ON], med[OFF], center=False)["on_off"],
"no_median": olfdyn.onoff_angles(r[ON], r[OFF], center=True)["on_off"],
"uncentred_no_median": olfdyn.onoff_angles(r[ON], r[OFF], center=False)["on_off"],
"zscored": olfdyn.onoff_angles(olfdyn.median_subtracted(rz)[ON], olfdyn.median_subtracted(rz)[OFF])["on_off"],
"zscored_pre": olfdyn.onoff_angles(olfdyn.median_subtracted(resp_zpre[s["name"], odor])[ON],
olfdyn.median_subtracted(resp_zpre[s["name"], odor])[OFF])["on_off"],
"frames": olfdyn.onoff_angles(olfdyn.median_subtracted(resp_fr[s["name"], odor])[ON],
olfdyn.median_subtracted(resp_fr[s["name"], odor])[OFF])["on_off"],
# centring controls: the same centred procedure between two windows of one period
"ctrl_on_halves": olfdyn.onoff_angles(med[32:72], med[72:112], center=True)["on_off"],
"ctrl_off_halves": olfdyn.onoff_angles(med[112:152], med[152:192], center=True)["on_off"],
"ctrl_pre_on": olfdyn.onoff_angles(med[PRE], med[ON], center=True)["on_off"],
"random": olfdyn.onoff_angles(noise[:80], noise[80:])["on_off"]})
ang = pd.DataFrame(rows)
print(ang.drop(columns=["ctrl_on_halves", "ctrl_off_halves", "ctrl_pre_on"]).round(1).to_string(index=False))
cols = ["on_off", "on_on", "off_off", "uncentred", "no_median", "uncentred_no_median", "zscored", "zscored_pre",
"frames", "ctrl_on_halves", "ctrl_off_halves", "ctrl_pre_on", "random"]
summary = pd.DataFrame({c: [ang[c].mean(), ang[c].std(ddof=1) / np.sqrt(len(ang)), ang[c].min(), ang[c].max(),
int((ang[c] > 90).sum())] for c in cols},
index=["mean", "sem_28", "min", "max", "n_above_90"]).T
print("\nOver the 28 session-odor values:\n", summary.round(2).to_string())
print("\nMean ON-OFF angle per odor (over sessions):", ang.groupby("odor").on_off.mean().round(1).to_dict())
print("Mean ON-OFF angle per mouse:", ang.groupby("mouse").on_off.mean().round(1).to_dict())
print("Mean ON-OFF angle per session:", ang.groupby("session").on_off.mean().round(1).to_dict())
print("SEM over the 4 odor means:", round(ang.groupby("odor").on_off.mean().std(ddof=1) / 2, 2))
cos = []
for (name, odor), r in resp.items():
med = olfdyn.median_subtracted(r)
x = np.vstack([med[ON], med[OFF]]); x = x - x.mean(axis=0)
a, b = x[:80].mean(0), x[80:].mean(0)
cos.append(a @ b / np.linalg.norm(a) / np.linalg.norm(b))
print("Cosine between the mean ON and mean OFF vectors after centring over the 160 frames: min", round(min(cos), 6), "max", round(max(cos), 6))
print("Expected 121.8 +/- 2.98")
session mouse odor on_off on_on off_off uncentred no_median uncentred_no_median zscored zscored_pre frames random
sub-134 ses-22 134 msc 138.2 42.7 41.6 78.4 126.8 51.0 129.4 135.7 129.5 90.4
sub-134 ses-22 134 ace 153.9 23.6 29.6 52.4 145.4 52.6 138.1 152.1 150.3 90.3
sub-134 ses-22 134 als 110.0 72.3 66.2 80.1 135.1 90.9 101.0 101.1 104.2 90.2
sub-134 ses-22 134 eug 143.1 33.9 41.4 81.2 131.6 72.0 134.1 145.7 138.8 90.3
sub-141 ses-1 141 msc 143.4 41.3 33.1 58.9 136.5 45.7 132.2 146.7 145.8 90.3
sub-141 ses-1 141 ace 146.6 35.1 29.9 45.5 140.1 44.3 132.8 145.8 147.6 90.4
sub-141 ses-1 141 als 102.5 79.0 81.1 30.0 101.7 18.5 102.7 103.5 103.3 90.4
sub-141 ses-1 141 eug 146.1 37.3 31.0 45.5 140.1 44.2 127.7 140.0 146.3 90.4
sub-141 ses-5 141 msc 143.7 36.0 38.0 78.9 136.6 65.5 133.9 143.8 143.3 90.4
sub-141 ses-5 141 ace 153.5 24.3 28.8 59.7 146.5 55.4 141.0 155.0 153.4 90.4
sub-141 ses-5 141 als 128.9 49.0 54.4 75.1 129.5 63.2 119.6 127.8 128.2 90.3
sub-141 ses-5 141 eug 134.9 42.8 49.6 81.7 128.6 72.4 124.0 131.2 134.4 90.3
sub-141 ses-8 141 msc 144.2 37.9 34.4 69.2 132.6 59.9 127.9 149.6 144.6 90.3
sub-141 ses-8 141 ace 135.2 42.9 45.8 71.0 126.2 44.8 124.4 133.7 135.6 90.4
sub-141 ses-8 141 als 127.2 55.0 50.9 59.7 123.5 26.5 126.3 124.0 125.6 90.4
sub-141 ses-8 141 eug 130.0 52.4 47.6 79.6 128.4 60.9 118.1 120.9 128.6 90.4
sub-148 ses-2 148 msc 112.0 70.6 66.8 77.8 107.3 74.0 108.2 112.5 109.5 90.3
sub-148 ses-2 148 ace 116.2 71.5 58.0 77.4 113.1 67.4 113.4 116.9 113.5 90.4
sub-148 ses-2 148 als 102.8 81.9 75.5 77.5 99.6 64.7 101.5 102.6 100.9 90.4
sub-148 ses-2 148 eug 110.1 71.7 72.1 81.8 108.2 83.3 106.3 107.1 107.0 90.3
sub-148 ses-3 148 msc 112.3 71.6 68.3 79.9 108.9 81.2 112.3 114.2 110.1 90.3
sub-148 ses-3 148 ace 138.5 43.8 38.1 77.6 138.8 80.5 129.9 136.5 137.5 90.4
sub-148 ses-3 148 als 101.7 85.8 76.9 51.8 101.0 56.0 101.9 102.4 100.6 90.3
sub-148 ses-3 148 eug 134.0 48.6 46.1 89.6 128.5 88.6 121.9 126.0 131.5 90.3
sub-149 ses-1 149 msc 105.6 77.2 77.4 84.5 107.6 68.1 106.2 107.3 104.2 90.3
sub-149 ses-1 149 ace 136.4 46.6 41.1 92.7 132.7 73.4 115.6 121.2 133.9 90.3
sub-149 ses-1 149 als 104.0 78.5 76.4 82.1 101.6 68.1 105.1 105.3 102.5 90.4
sub-149 ses-1 149 eug 125.5 58.5 53.2 90.7 122.0 86.3 110.0 113.1 120.9 90.4
Over the 28 session-odor values:
mean sem_28 min max n_above_90
on_off 127.88 3.20 101.69 153.91 28.0
on_on 53.99 3.51 23.56 85.79 0.0
off_off 51.90 3.22 28.78 81.12 0.0
uncentred 71.80 2.93 30.04 92.74 2.0
no_median 124.24 2.75 99.63 146.48 28.0
uncentred_no_median 62.84 3.36 18.49 90.93 1.0
zscored 119.48 2.34 101.01 140.98 28.0
zscored_pre 125.78 3.27 101.14 154.96 28.0
frames 126.11 3.29 100.63 153.40 28.0
ctrl_on_halves 101.13 0.98 91.99 111.62 28.0
ctrl_off_halves 101.63 0.72 96.36 113.46 28.0
ctrl_pre_on 129.72 3.26 101.65 157.43 28.0
random 90.35 0.01 90.25 90.44 28.0
Mean ON-OFF angle per odor (over sessions): {'ace': 140.0, 'als': 111.0, 'eug': 131.9, 'msc': 128.5}
Mean ON-OFF angle per mouse: {'134': 136.3, '141': 136.3, '148': 116.0, '149': 117.9}
Mean ON-OFF angle per session: {'sub-134 ses-22': 136.3, 'sub-141 ses-1': 134.7, 'sub-141 ses-5': 140.3, 'sub-141 ses-8': 134.1, 'sub-148 ses-2': 110.3, 'sub-148 ses-3': 121.6, 'sub-149 ses-1': 117.9}
SEM over the 4 odor means: 6.12
Cosine between the mean ON and mean OFF vectors after centring over the 160 frames: min -1.0 max -1.0
Expected 121.8 +/- 2.98
fig, ax = plt.subplots(figsize=(7, 3))
for i, c in enumerate(["on_off", "uncentred", "no_median", "zscored", "random"]):
ax.scatter(np.full(len(ang), i) + np.linspace(-0.2, 0.2, len(ang)), ang[c], s=12)
ax.axhline(90, color="grey", lw=0.8)
ax.axhline(121.8, color="k", lw=0.8, ls="--")
ax.set_xticks(range(5), ["main", "uncentred", "no median", "z-scored", "random"])
ax.set_ylabel("median ON-OFF angle (deg)")
ax.set_ylim(0, 180)
plt.tight_layout()
print("All plotted values lie in 0-180 degrees:", bool(ang[["on_off", "uncentred", "no_median", "zscored", "random"]].stack().between(0, 180).all()))
All plotted values lie in 0-180 degrees: True
# Synthetic check of the centring: ON and OFF frames drawn around two fixed patterns that are
# close in angle (same mean pattern plus a small difference), with frame noise.
rng_syn = np.random.default_rng(21)
rows = []
for sep in [0.1, 0.2, 0.5]:
base_v = rng_syn.standard_normal(200) + 2.0
diff = rng_syn.standard_normal(200) * sep
on_s = base_v + diff / 2 + 0.3 * rng_syn.standard_normal((80, 200))
off_s = base_v - diff / 2 + 0.3 * rng_syn.standard_normal((80, 200))
rows.append({"pattern_angle": float(olfdyn.angles((base_v + diff / 2)[None], (base_v - diff / 2)[None])[0, 0]),
"uncentred": olfdyn.onoff_angles(on_s, off_s, center=False)["on_off"],
"centred": olfdyn.onoff_angles(on_s, off_s, center=True)["on_off"]})
print(pd.DataFrame(rows).round(1).to_string(index=False))
pattern_angle uncentred centred
2.4 11.3 91.6
4.9 11.6 96.1
13.0 16.6 115.6
Result: partially reproduced. With the code's procedure the median ON-OFF angle exceeds 90 degrees in all 28 session-odor pairs (101.7 to 153.9) and in each odor's mean (111.0 to 140.0), but the mean, 127.9 degrees (SEM 3.2 over the 28 values), is outside 121.8 +/- 2.98. Without the median subtraction the mean is 124.2, and with per-trial z-scored traces 119.5, both within tolerance; with F z-scored on the 2 s before onset it is 125.8, and with windows counted in frames as the code does 126.1, both outside. Random vectors give 90.4, and the ON-ON and OFF-OFF medians average 54.0 and 51.9 degrees.
The angle above 90 degrees does not show that ON and OFF are inverted. After centring, the mean ON vector is exactly opposite the mean OFF vector (cosine -1 in every case), and without centring the mean angle is 71.8 degrees, with 2 of 28 above 90. The same centred procedure gives 101.1 degrees between the first and second halves of the ON window (28 of 28 above 90), 101.6 between the halves of the OFF window and 129.7 between the 2 s before onset and the ON window. On synthetic frames around two patterns 2.4, 4.9 and 13.0 degrees apart, it gives 91.6, 96.1 and 115.6 degrees, while the uncentred angles are 11.3, 11.6 and 16.6.
fig2f-onoff-anticorrelated-mouse¶
Time-averaged glomerular ON and OFF responses to the same odor are anticorrelated; glomeruli active during the odor tend to be suppressed after it and the reverse.
Result: not reproduced. The time-averaged ON and OFF vectors are positively correlated in all 28 session-odor pairs (median r = 0.481, Wilcoxon p = 7.5e-09), for every odor (medians 0.391 to 0.535) and every mouse (0.245 to 0.596). The values without the median subtraction are identical by construction, as the subtraction shifts every glomerulus of a vector by the same amount, so they are not a separate check. With F z-scored on the 2 s before onset the median r is 0.526 (0 of 28 below 0), and with windows counted in frames 0.446 (0 of 28). With per-trial z-scored traces, the normalisation closest to the paper's deposited arrays, the median r is 0.068 (12 of 28 below 0, p = 0.33). That normalisation can produce negative correlations by itself: on synthetic glomeruli with excitation only, z-scoring each trial over its window gives r = -0.90 between ON and 10 to 15 s after offset and -0.93 at 45 to 50 s, against 0.30 and 0.13 with the pre-onset z-score.
fig4b-reference-odor-mouse¶
In the mouse OB the ON response to a reference odor is closer in angle to the ON responses of the other odors than to any OFF response, including its own, and this holds for every odor of the panel.
def reference_angles(name, top=3):
blocks = []
for odor in olfdyn.BLOCK_ODORS:
med = olfdyn.median_subtracted(resp[name, odor])
blocks += [med[ON], med[OFF]]
x = np.vstack(blocks)
x = x - x.mean(axis=0)
if top:
u, sv, vt = np.linalg.svd(x, full_matrices=False)
x = x @ vt[:top].T
means = np.stack([x[i * 80:(i + 1) * 80].mean(0) for i in range(8)]) # odor1 ON, odor1 OFF, ...
on, off = means[0::2], means[1::2]
out = []
for i, odor in enumerate(olfdyn.BLOCK_ODORS):
a_on = olfdyn.angles(on[i:i + 1], on)[0]
a_off = olfdyn.angles(on[i:i + 1], off)[0]
others = np.delete(a_on, i)
out.append({"session": name, "ref": odor, "on_on_mean": others.mean(), "on_off_mean": a_off.mean(),
"on_on_max": others.max(), "on_off_min": a_off.min(), "own_off": a_off[i]})
return out
for top, label in [(3, "top three PCs (main)"), (None, "full dimension (variant)")]:
ra = pd.DataFrame([r for s in block for r in reference_angles(s["name"], top)])
print(f"--- {label}")
print(ra.round(1).to_string(index=False))
for odor, g in ra.groupby("ref", sort=False):
t = stats.ttest_rel(g.on_on_mean, g.on_off_mean)
print(f"ref {odor}: ON-ON {g.on_on_mean.mean():.1f}, ON-OFF {g.on_off_mean.mean():.1f}, paired t p = {t.pvalue:.2g}; "
f"sessions where max ON-ON < min ON-OFF: {int((g.on_on_max < g.on_off_min).sum())} of {len(g)}")
print(f"Session-mean ON-ON angles range from {ra.on_on_mean.min():.1f} to {ra.on_on_mean.max():.1f} degrees")
if top == 3:
ra_main = ra
--- top three PCs (main)
session ref on_on_mean on_off_mean on_on_max on_off_min own_off
sub-134 ses-22 msc 98.7 105.3 104.0 79.0 79.0
sub-134 ses-22 ace 117.7 129.6 151.5 62.9 62.9
sub-134 ses-22 als 115.3 37.2 151.5 13.1 13.1
sub-134 ses-22 eug 100.6 99.2 110.3 73.8 73.8
sub-141 ses-1 msc 98.3 118.7 118.5 88.0 88.0
sub-141 ses-1 ace 104.2 112.0 140.0 56.2 56.2
sub-141 ses-1 als 123.1 49.3 140.0 4.0 4.0
sub-141 ses-1 eug 106.6 99.5 110.6 49.3 49.3
sub-141 ses-5 msc 85.6 120.9 111.5 93.5 93.5
sub-141 ses-5 ace 113.5 119.7 119.4 53.5 53.5
sub-141 ses-5 als 107.8 78.0 127.4 62.2 62.2
sub-141 ses-5 eug 101.9 93.2 127.4 79.2 87.9
sub-141 ses-8 msc 113.2 122.0 140.8 79.1 79.1
sub-141 ses-8 ace 106.4 127.1 131.9 90.3 90.3
sub-141 ses-8 als 110.2 47.7 140.8 21.3 47.0
sub-141 ses-8 eug 90.6 92.8 112.6 79.2 79.2
sub-148 ses-2 msc 106.2 109.3 154.9 76.1 76.1
sub-148 ses-2 ace 105.5 126.1 144.8 94.2 94.2
sub-148 ses-2 als 106.1 71.9 144.8 35.2 35.2
sub-148 ses-2 eug 112.1 82.4 154.9 63.2 63.2
sub-148 ses-3 msc 116.3 62.2 122.6 23.1 23.1
sub-148 ses-3 ace 103.4 130.3 115.0 106.4 132.9
sub-148 ses-3 als 111.8 65.2 114.3 12.6 12.6
sub-148 ses-3 eug 104.4 115.0 122.6 99.4 114.4
sub-149 ses-1 msc 114.4 64.9 137.6 35.0 35.0
sub-149 ses-1 ace 89.2 147.5 132.1 126.0 126.0
sub-149 ses-1 als 110.1 68.2 131.9 35.1 35.1
sub-149 ses-1 eug 93.3 142.0 137.6 115.5 144.5
ref msc: ON-ON 104.7, ON-OFF 100.5, paired t p = 0.76; sessions where max ON-ON < min ON-OFF: 0 of 7
ref ace: ON-ON 105.7, ON-OFF 127.5, paired t p = 0.018; sessions where max ON-ON < min ON-OFF: 0 of 7
ref als: ON-ON 112.1, ON-OFF 59.6, paired t p = 0.00035; sessions where max ON-ON < min ON-OFF: 0 of 7
ref eug: ON-ON 101.4, ON-OFF 103.4, paired t p = 0.83; sessions where max ON-ON < min ON-OFF: 0 of 7
Session-mean ON-ON angles range from 85.6 to 123.1 degrees
--- full dimension (variant)
session ref on_on_mean on_off_mean on_on_max on_off_min own_off
sub-134 ses-22 msc 97.0 104.0 101.0 81.5 81.5
sub-134 ses-22 ace 111.6 122.7 133.1 66.3 66.3
sub-134 ses-22 als 107.9 66.5 133.1 56.6 56.6
sub-134 ses-22 eug 101.1 97.9 110.3 74.9 74.9
sub-141 ses-1 msc 100.2 108.9 117.0 92.3 92.3
sub-141 ses-1 ace 101.3 109.7 127.2 74.7 74.7
sub-141 ses-1 als 117.1 67.6 127.2 33.9 33.9
sub-141 ses-1 eug 105.0 98.6 107.5 61.4 61.4
sub-141 ses-5 msc 90.3 106.2 110.6 88.4 88.4
sub-141 ses-5 ace 111.8 113.3 116.8 70.2 70.2
sub-141 ses-5 als 104.6 84.4 116.8 72.3 72.3
sub-141 ses-5 eug 99.2 95.0 113.0 87.6 87.6
sub-141 ses-8 msc 103.8 116.1 114.1 84.8 84.8
sub-141 ses-8 ace 99.3 115.1 110.4 91.0 91.0
sub-141 ses-8 als 102.7 86.0 114.1 70.1 70.1
sub-141 ses-8 eug 98.7 91.5 111.1 79.5 79.5
sub-148 ses-2 msc 100.0 105.0 113.5 82.2 82.2
sub-148 ses-2 ace 98.4 110.3 110.9 92.6 92.6
sub-148 ses-2 als 102.9 88.7 110.9 69.5 69.5
sub-148 ses-2 eug 103.0 94.7 113.5 76.3 76.3
sub-148 ses-3 msc 107.6 83.6 112.0 59.9 59.9
sub-148 ses-3 ace 100.4 117.1 110.7 104.9 107.9
sub-148 ses-3 als 106.5 81.3 110.7 46.5 46.5
sub-148 ses-3 eug 100.3 105.1 112.0 100.2 108.4
sub-149 ses-1 msc 111.6 73.0 130.0 55.5 55.5
sub-149 ses-1 ace 96.7 130.4 129.9 117.2 117.2
sub-149 ses-1 als 107.6 77.7 126.3 63.0 63.0
sub-149 ses-1 eug 98.3 122.6 130.0 117.0 117.0
ref msc: ON-ON 101.5, ON-OFF 99.5, paired t p = 0.81; sessions where max ON-ON < min ON-OFF: 0 of 7
ref ace: ON-ON 102.8, ON-OFF 116.9, paired t p = 0.0096; sessions where max ON-ON < min ON-OFF: 0 of 7
ref als: ON-ON 107.0, ON-OFF 78.9, paired t p = 0.0013; sessions where max ON-ON < min ON-OFF: 0 of 7
ref eug: ON-ON 100.8, ON-OFF 100.8, paired t p = 0.99; sessions where max ON-ON < min ON-OFF: 0 of 7
Session-mean ON-ON angles range from 90.3 to 117.1 degrees
Result: partially reproduced. Reading (a) holds for one reference odor: the acetophenone ON vector is closer to the other ON vectors (105.7 degrees) than to the OFF vectors (127.5; paired p = 0.018). For methyl salicylate and eugenol the two do not differ (p = 0.76 and 0.83), and the allyl sulfide ON vector is closer to the OFF vectors (ON-ON 112.1, ON-OFF 59.6, p = 0.00035), the reverse of the claim, so it does not count for reading (a), which the t-test alone (two-sided) would not exclude. Reading (b) holds in none of the seven sessions for any reference odor. The full-dimensional variant gives the same pattern (acetophenone p = 0.0096 in the claimed direction, allyl sulfide p = 0.0013 reversed). The session-mean ON-ON angles here are 85.6 to 123.1 degrees; the published version of the paper gives a median of 44.7 degrees for the acetophenone reference (iScience, Figure 2G), so its procedure probably differs from this one.
fig4d-clustering-mouse¶
Hierarchical clustering of glomerular responses cleanly separates ON responses from OFF responses.
rows = []
for s in block:
for label, src in [("dff", resp), ("zscored", resp_z)]:
vec = np.vstack([np.vstack([src[s["name"], o][ON].mean(0), src[s["name"], o][OFF].mean(0)])
for o in olfdyn.BLOCK_ODORS]) # ON, OFF per odor
cl = fcluster(linkage(pdist(vec, "cosine"), "complete"), 2, "maxclust")
on_cl, off_cl = set(cl[0::2]), set(cl[1::2])
rows.append({"session": s["name"], "mouse": s["subject"], "rois": s["n_rois"], "vectors": label,
"clusters_on": sorted(on_cl), "clusters_off": sorted(off_cl),
"separates": len(on_cl) == 1 and len(off_cl) == 1 and on_cl != off_cl})
cl = pd.DataFrame(rows)
print(cl.to_string(index=False))
for label, g in cl.groupby("vectors"):
print(f"{label}: ON and OFF separated in {int(g.separates.sum())} of {len(g)} sessions; "
f"198-ROI session: {bool(g[g.rois == 198].separates.iloc[0])}")
session mouse rois vectors clusters_on clusters_off separates sub-134 ses-22 134 73 dff [1, 2] [1] False sub-134 ses-22 134 73 zscored [1, 2] [1] False sub-141 ses-1 141 210 dff [1, 2] [1, 2] False sub-141 ses-1 141 210 zscored [1, 2] [1] False sub-141 ses-5 141 153 dff [1, 2] [1] False sub-141 ses-5 141 153 zscored [1, 2] [1] False sub-141 ses-8 141 198 dff [1, 2] [1] False sub-141 ses-8 141 198 zscored [1, 2] [1] False sub-148 ses-2 148 286 dff [1, 2] [1, 2] False sub-148 ses-2 148 286 zscored [1, 2] [1] False sub-148 ses-3 148 318 dff [1, 2] [1, 2] False sub-148 ses-3 148 318 zscored [1, 2] [2] False sub-149 ses-1 149 219 dff [1, 2] [2] False sub-149 ses-1 149 219 zscored [1, 2] [1] False dff: ON and OFF separated in 0 of 7 sessions; 198-ROI session: False zscored: ON and OFF separated in 0 of 7 sessions; 198-ROI session: False
Result: not reproduced. In none of the seven sessions, the 198-glomerulus session included, does the two-cluster cut of the dendrogram put the four ON vectors in one cluster and the four OFF vectors in the other, with dF/F or with z-scored traces (0 of 7 each).
fig6e-contrast-msc¶
Methyl salicylate presented 1 s after another odor evokes a response less correlated with that odor's solitary response than the solitary methyl salicylate response is (history dependent contrast enhancement; one-tailed t-test p < 0.01, 3 sessions from 2 mice).
SEQ_B = slice(128, 208) # the second odor's 5 s, 6 to 11 s after the pair's onset
rows = []
for s in block:
seq = olfdyn.dff_pre(s["seq"])
seq_z = olfdyn.zscore_trials(s["seq"][:, :912])
seq_zp = olfdyn.zscore_pre(s["seq"])
lag = (s["pres"][s["pres"].kind == "seq_second"].onset.to_numpy()
- s["pres"][s["pres"].kind == "seq_first"].onset.to_numpy())
msc_sol = resp[s["name"], "msc"][ON].mean(0)
for x, lag_x in zip(s["seq_pres"].odor.unique(), [lag[(s["seq_pres"].odor == o).to_numpy()].mean() for o in s["seq_pres"].odor.unique()]):
sel = (s["seq_pres"].odor == x).to_numpy()
x_sol = resp[s["name"], x][ON].mean(0)
msc_seq = seq[sel].mean(axis=0)[SEQ_B].mean(0)
x_seq_first = seq[sel].mean(axis=0)[ON].mean(0)
zx = resp_z[s["name"], x][ON].mean(0)
zm = resp_z[s["name"], "msc"][ON].mean(0)
zs = seq_z[sel].mean(axis=0)[SEQ_B].mean(0)
zpx = resp_zpre[s["name"], x][ON].mean(0)
zpm = resp_zpre[s["name"], "msc"][ON].mean(0)
zps = seq_zp[sel].mean(axis=0)[SEQ_B].mean(0)
rows.append({"session": s["name"], "mouse": s["subject"], "rois": s["n_rois"], "first": x, "pairs": int(sel.sum()), "lag_s": lag_x,
"r_solitary": olfdyn.corr(x_sol, msc_sol), "r_sequential": olfdyn.corr(x_sol, msc_seq),
"r_msc_sol_vs_seq": olfdyn.corr(msc_sol, msc_seq),
"r_first_sol_vs_seq": olfdyn.corr(x_sol, x_seq_first),
"z_solitary": olfdyn.corr(zx, zm), "z_sequential": olfdyn.corr(zx, zs),
"zpre_solitary": olfdyn.corr(zpx, zpm), "zpre_sequential": olfdyn.corr(zpx, zps)})
ce = pd.DataFrame(rows)
print(ce.round(3).to_string(index=False))
d = ce.r_solitary - ce.r_sequential
t = stats.ttest_rel(ce.r_solitary, ce.r_sequential, alternative="greater")
print(f"\nAll 21 session-odor pairs: solitary r {ce.r_solitary.mean():.3f}, sequential r {ce.r_sequential.mean():.3f}, "
f"difference {d.mean():.3f} (positive in {int((d > 0).sum())} of 21), one-tailed paired t p = {t.pvalue:.2g}")
ps = ce.groupby("session")[["r_solitary", "r_sequential"]].mean()
t2 = stats.ttest_rel(ps.r_solitary, ps.r_sequential, alternative="greater")
print(f"Per session (mean over first odors, n=7): difference {(ps.r_solitary - ps.r_sequential).mean():.3f}, p = {t2.pvalue:.2g}")
for x, g in ce.groupby("first"):
tx = stats.ttest_rel(g.r_solitary, g.r_sequential, alternative="greater")
print(f"first odor {x}: difference {(g.r_solitary - g.r_sequential).mean():.3f}, p = {tx.pvalue:.2g}")
tz = stats.ttest_rel(ce.z_solitary, ce.z_sequential, alternative="greater")
print(f"Variant, per-trial z-scored traces: solitary {ce.z_solitary.mean():.3f}, sequential {ce.z_sequential.mean():.3f}, "
f"one-tailed paired p = {tz.pvalue:.2g}")
tzp = stats.ttest_rel(ce.zpre_solitary, ce.zpre_sequential, alternative="greater")
print(f"Variant, z-scored with the mean and SD of the 2 s before the pair's onset: solitary {ce.zpre_solitary.mean():.3f}, "
f"sequential {ce.zpre_sequential.mean():.3f}, one-tailed paired p = {tzp.pvalue:.2g}")
# The published version's Table 1 (iScience 2026) names the sessions of this panel by glomerulus count: 286, 219, 318.
t3 = ce[ce.rois.isin([286, 219, 318])]
tt3 = stats.ttest_rel(t3.r_solitary, t3.r_sequential, alternative="greater")
tz3 = stats.ttest_rel(t3.z_solitary, t3.z_sequential, alternative="greater")
print(f"Variant, the three sessions with 286, 219 and 318 ROIs ({sorted(t3.session.unique())}, {len(t3)} pairs): "
f"solitary {t3.r_solitary.mean():.3f}, sequential {t3.r_sequential.mean():.3f}, one-tailed p = {tt3.pvalue:.2g}; "
f"per-trial z-scored p = {tz3.pvalue:.2g}")
pm = ce.groupby("mouse")[["r_solitary", "r_sequential"]].mean()
print("Per mouse (mean):\n", pm.round(3).to_string())
session mouse rois first pairs lag_s r_solitary r_sequential r_msc_sol_vs_seq r_first_sol_vs_seq z_solitary z_sequential zpre_solitary zpre_sequential
sub-134 ses-22 134 73 ace 10 6.046 0.390 0.654 0.834 0.985 0.606 0.527 0.333 0.752
sub-134 ses-22 134 73 als 10 6.043 0.135 0.214 0.945 0.772 0.341 0.311 0.266 0.398
sub-134 ses-22 134 73 eug 10 6.043 0.082 0.121 0.924 0.885 0.169 0.002 0.064 0.052
sub-141 ses-1 141 210 ace 10 6.041 0.628 0.652 0.923 0.981 0.533 0.415 0.596 0.626
sub-141 ses-1 141 210 als 10 6.041 0.462 0.177 0.898 0.210 0.379 0.132 0.358 0.166
sub-141 ses-1 141 210 eug 10 6.040 0.321 0.329 0.858 0.973 0.367 0.214 0.278 0.277
sub-141 ses-5 141 153 ace 10 6.047 0.062 0.249 0.838 0.943 0.199 0.023 0.088 0.314
sub-141 ses-5 141 153 als 10 6.044 0.394 0.375 0.760 0.903 0.360 0.079 0.305 0.427
sub-141 ses-5 141 153 eug 10 6.042 0.497 0.334 0.812 0.903 0.463 0.278 0.617 0.465
sub-141 ses-8 141 198 ace 10 6.047 0.365 0.377 0.858 0.942 0.467 0.228 0.407 0.459
sub-141 ses-8 141 198 als 10 6.042 0.143 0.338 0.842 0.755 0.088 0.000 0.097 0.230
sub-141 ses-8 141 198 eug 10 6.041 0.131 0.215 0.926 0.897 0.210 0.134 0.085 0.149
sub-148 ses-2 148 286 ace 10 6.045 0.317 0.126 0.650 0.573 0.422 0.098 0.360 0.125
sub-148 ses-2 148 286 als 10 6.042 0.390 0.357 0.820 0.454 0.196 0.123 0.352 0.308
sub-148 ses-2 148 286 eug 10 6.048 0.018 -0.022 0.863 0.760 0.175 0.003 0.049 0.056
sub-148 ses-3 148 318 ace 10 6.045 0.413 0.170 0.544 0.916 0.564 0.145 0.493 0.302
sub-148 ses-3 148 318 als 10 6.044 0.128 0.325 0.677 0.832 0.134 0.137 0.229 0.328
sub-148 ses-3 148 318 eug 10 6.045 0.226 0.125 0.721 0.839 0.392 0.048 0.275 0.157
sub-149 ses-1 149 219 ace 10 6.043 -0.193 -0.343 0.633 0.938 -0.009 -0.161 -0.051 -0.049
sub-149 ses-1 149 219 als 10 6.043 0.116 -0.004 0.660 0.573 0.151 -0.111 0.204 -0.036
sub-149 ses-1 149 219 eug 10 6.043 -0.220 -0.159 0.758 0.885 -0.132 -0.140 -0.136 -0.084
All 21 session-odor pairs: solitary r 0.229, sequential r 0.219, difference 0.009 (positive in 10 of 21), one-tailed paired t p = 0.39
Per session (mean over first odors, n=7): difference 0.009, p = 0.4
first odor ace: difference 0.014, p = 0.43
first odor als: difference -0.002, p = 0.51
first odor eug: difference 0.016, p = 0.33
Variant, per-trial z-scored traces: solitary 0.289, sequential 0.118, one-tailed paired p = 5.6e-07
Variant, z-scored with the mean and SD of the 2 s before the pair's onset: solitary 0.251, sequential 0.258, one-tailed paired p = 0.58
Variant, the three sessions with 286, 219 and 318 ROIs (['sub-148 ses-2', 'sub-148 ses-3', 'sub-149 ses-1'], 9 pairs): solitary 0.133, sequential 0.064, one-tailed p = 0.082; per-trial z-scored p = 0.0025
Per mouse (mean):
r_solitary r_sequential
mouse
134 0.202 0.330
141 0.334 0.338
148 0.249 0.180
149 -0.099 -0.169
Result: not reproduced. Over the 21 session-odor pairs the correlation of solitary X with solitary methyl salicylate (mean 0.229) does not exceed its correlation with methyl salicylate in the (X)msc pairs (0.219; difference 0.009, positive in 10 of 21, one-tailed p = 0.39). Per session p = 0.40, and per first odor p = 0.33 to 0.51. The sequential methyl salicylate response stays close to the solitary one (r = 0.544 to 0.945). On the three sessions the published version's Table 1 names for this panel (286, 219 and 318 ROIs; sub-148 ses-2, sub-148 ses-3, sub-149 ses-1) the solitary r is 0.133 and the sequential r 0.064 (p = 0.082). With F z-scored on the 2 s before the pair's onset there is no difference (0.251 against 0.258, p = 0.58). With per-trial z-scored traces the solitary r is 0.289 and the sequential r 0.118 (p = 5.6e-07; p = 0.0025 on the three sessions), as the claim states; the synthetic check under fig2f shows that this normalisation alone can turn a response into negative correlations with later activity.
fig7d-persistent-anticorrelation¶
After odor offset the glomerular population activity stays negatively correlated with the odor-evoked ON response for several tens of seconds, although overall activity returns to baseline quickly.
BINS = [slice(112 + 80 * i, 192 + 80 * i) for i in range(10)] # 0-5 s ... 45-50 s after offset
rows = []
for s in block:
for odor in olfdyn.BLOCK_ODORS:
r = resp[s["name"], odor]
med = olfdyn.median_subtracted(r)
sel = (s["sol_pres"].odor == odor).to_numpy()
f0 = s["pre_block"][odor][0]
pb = ((s["sol"][sel] - f0) / f0).mean(axis=0).astype(np.float64)
pbm = olfdyn.median_subtracted(pb)
rz = olfdyn.median_subtracted(resp_z[s["name"], odor])
rzp = olfdyn.median_subtracted(resp_zpre[s["name"], odor])
rfr = olfdyn.median_subtracted(resp_fr[s["name"], odor])
for i, b in enumerate(BINS):
rows.append({"session": s["name"], "mouse": s["subject"], "odor": odor, "bin_s": f"{5 * i}-{5 * i + 5}",
"r_main": olfdyn.corr(med[ON].mean(0), med[b].mean(0)),
"r_no_median": olfdyn.corr(r[ON].mean(0), r[b].mean(0)),
"r_preblock_f0": olfdyn.corr(pbm[ON].mean(0), pbm[b].mean(0)),
"r_zscored": olfdyn.corr(rz[ON].mean(0), rz[b].mean(0)),
"r_zscored_pre": olfdyn.corr(rzp[ON].mean(0), rzp[b].mean(0)),
"r_frames": olfdyn.corr(rfr[ON].mean(0), rfr[b].mean(0)),
"mean_dff": float(r[b].mean()), "mean_dff_on": float(r[ON].mean())})
pc = pd.DataFrame(rows)
order = [f"{5 * i}-{5 * i + 5}" for i in range(10)]
out = []
for b in order:
g = pc[pc.bin_s == b]
row = {"bin_s": b}
for c in ["r_main", "r_preblock_f0", "r_zscored", "r_zscored_pre", "r_frames"]:
row[c + "_median"] = g[c].median()
row[c + "_p"] = stats.wilcoxon(g[c]).pvalue
row[c + "_neg"] = int((g[c] < 0).sum())
row["mean_dff"] = g.mean_dff.mean()
out.append(row)
print("r_no_median equals r_main by construction:", bool(np.allclose(pc.r_main, pc.r_no_median)))
print("Median over 28 session-odor values per bin after offset, two-sided Wilcoxon p, and how many are below 0:")
print(pd.DataFrame(out).round(4).to_string(index=False))
print("\nMean dF/F during ON (over sessions and odors):", round(pc.mean_dff_on.mean(), 4))
print("\nr_main median per mouse and bin:")
print(pc.pivot_table(index="mouse", columns="bin_s", values="r_main", aggfunc="median")[order].round(3).to_string())
print("\nPre-block baseline length (s) per session and odor:",
{s["name"]: {o: round(v[1]) for o, v in s["pre_block"].items()} for s in block})
r_no_median equals r_main by construction: True
Median over 28 session-odor values per bin after offset, two-sided Wilcoxon p, and how many are below 0:
bin_s r_main_median r_main_p r_main_neg r_preblock_f0_median r_preblock_f0_p r_preblock_f0_neg r_zscored_median r_zscored_p r_zscored_neg r_zscored_pre_median r_zscored_pre_p r_zscored_pre_neg r_frames_median r_frames_p r_frames_neg mean_dff
0-5 0.4810 0.0000 0 0.5998 0.0 0 0.0679 0.3272 12 0.5265 0.0000 0 0.4461 0.0000 0 0.0707
5-10 0.3116 0.0000 2 0.4342 0.0 1 -0.3438 0.0000 28 0.3116 0.0000 0 0.3256 0.0000 2 0.0450
10-15 0.3283 0.0000 3 0.4276 0.0 1 -0.4124 0.0000 28 0.3127 0.0000 2 0.3385 0.0000 3 0.0265
15-20 0.3254 0.0000 2 0.4297 0.0 2 -0.4796 0.0000 28 0.2935 0.0000 1 0.3324 0.0000 3 0.0192
20-25 0.2524 0.0000 4 0.3932 0.0 2 -0.4712 0.0000 28 0.2634 0.0000 3 0.2321 0.0000 4 0.0111
25-30 0.1803 0.0004 7 0.3814 0.0 3 -0.4656 0.0000 28 0.2552 0.0000 3 0.2522 0.0001 5 0.0062
30-35 0.1509 0.0001 6 0.3856 0.0 3 -0.4574 0.0000 28 0.2099 0.0000 5 0.1724 0.0004 6 0.0049
35-40 0.1563 0.0041 8 0.3254 0.0 4 -0.4779 0.0000 27 0.1519 0.0012 5 0.1609 0.0014 7 0.0069
40-45 0.1795 0.0006 9 0.3461 0.0 2 -0.4529 0.0000 26 0.1752 0.0001 4 0.1515 0.0052 9 0.0068
45-50 0.0943 0.0136 10 0.3037 0.0 3 -0.4224 0.0000 27 0.1770 0.0005 8 0.2033 0.0003 7 0.0045
Mean dF/F during ON (over sessions and odors): 0.0539
r_main median per mouse and bin:
bin_s 0-5 5-10 10-15 15-20 20-25 25-30 30-35 35-40 40-45 45-50
mouse
134 0.508 0.194 0.158 0.100 0.080 0.040 0.053 0.034 -0.063 -0.019
141 0.596 0.359 0.199 0.185 0.176 0.032 0.070 0.007 0.047 -0.054
148 0.424 0.299 0.385 0.374 0.342 0.304 0.298 0.212 0.283 0.273
149 0.245 0.274 0.370 0.356 0.343 0.288 0.338 0.307 0.226 0.263
Pre-block baseline length (s) per session and odor: {'sub-134 ses-22': {'msc': 19, 'ace': 900, 'als': 900, 'eug': 900}, 'sub-141 ses-1': {'msc': 24, 'ace': 900, 'als': 900, 'eug': 900}, 'sub-141 ses-5': {'msc': 14, 'ace': 900, 'als': 900, 'eug': 900}, 'sub-141 ses-8': {'msc': 20, 'ace': 900, 'als': 900, 'eug': 900}, 'sub-148 ses-2': {'msc': 848, 'ace': 900, 'als': 900, 'eug': 900}, 'sub-148 ses-3': {'msc': 848, 'ace': 900, 'als': 900, 'eug': 900}, 'sub-149 ses-1': {'msc': 848, 'ace': 900, 'als': 900, 'eug': 900}}
Result: not reproduced. After offset the population vector stays positively correlated with the ON vector: the median r over the 28 session-odor pairs is 0.481 at 0 to 5 s, 0.312 to 0.328 from 5 to 20 s, and 0.094 at 45 to 50 s, above 0 in every bin (p <= 0.014). The mean dF/F falls from 0.071 in the first 5 s after offset to 0.005 at 45 to 50 s, so the activity after offset is a fading positive copy of the ON pattern; in mice 148 and 149 r stays near 0.2 to 0.4 for the whole 50 s. With F0 from the 15 min before each block r is positive in every bin (medians 0.304 to 0.600), with F z-scored on the 2 s before onset too (0.152 to 0.527), and with windows counted in frames too (0.152 to 0.446). With per-trial z-scored traces r is negative from 5 s on (medians -0.344 to -0.480, 26 to 28 of 28 below 0), as the claim states; on synthetic glomeruli with excitation only, that normalisation gives negative values of the same sign and larger size (-0.90 and -0.93, fig2f). The values without the median subtraction equal the main ones by construction.
fig7h-iti-knn-mouse¶
Glomerular activity averaged over the 50 s after each odor offset classifies the preceding odor with high accuracy (k-NN, k = 20).
ITI = slice(112, 912) # 0 to 50 s after offset
rng = np.random.default_rng(7)
rows, null_main = [], None
for s in block:
x = olfdyn.dff_session(s["sol"][:, ITI], s).mean(axis=1).astype(np.float64) # 40 x ROIs
y = s["sol_pres"].odor.to_numpy()
nb = olfdyn.knn_neighbours(x, np.arange(len(x)), 20, "cosine")
acc = (olfdyn.knn_predict(nb, y) == y).mean()
null = olfdyn.permuted_accuracy(nb, y, np.arange(len(x)), rng, 1000)
xc = x - x.mean(0)
pcs = xc @ np.linalg.svd(xc, full_matrices=False)[2][:3].T
nb3 = olfdyn.knn_neighbours(pcs, np.arange(len(x)), 20, "cosine")
acc3 = (olfdyn.knn_predict(nb3, y) == y).mean()
small = {}
for k in [5, 9]: # at or below the 9 other trials of the own odor
nbk = olfdyn.knn_neighbours(x, np.arange(len(x)), k, "cosine")
nullk = olfdyn.permuted_accuracy(nbk, y, np.arange(len(x)), rng, 1000)
small[f"acc_k{k}"] = (olfdyn.knn_predict(nbk, y) == y).mean()
small[f"null_mean_k{k}"] = nullk.mean()
small[f"null_p99_k{k}"] = np.percentile(nullk, 99)
rows.append({"session": s["name"], "rois": s["n_rois"], "acc_full_k20": acc, "null_p99": np.percentile(null, 99),
"null_mean": null.mean(), "acc_top3pc_k20": acc3, **small,
"per_odor": {o: round(float((olfdyn.knn_predict(nb, y) == y)[y == o].mean()), 2) for o in olfdyn.BLOCK_ODORS}})
ik = pd.DataFrame(rows)
print(ik.drop(columns=[c for c in ik.columns if c.endswith(("_k5", "_k9")) or c in ("acc_k5", "acc_k9")]).round(3).to_string(index=False))
print("\nVariant, k at or below the class size (9 other trials of the own odor):")
print(ik[["session", "rois", "acc_k5", "null_mean_k5", "null_p99_k5", "acc_k9", "null_mean_k9", "null_p99_k9"]].round(3).to_string(index=False))
print("Chance 0.25. Main session: the one with 318 ROIs.")
session rois acc_full_k20 null_p99 null_mean acc_top3pc_k20 per_odor
sub-134 ses-22 73 0.225 0.400 0.169 0.250 {'msc': 0.0, 'ace': 0.9, 'als': 0.0, 'eug': 0.0}
sub-141 ses-1 210 0.250 0.375 0.169 0.275 {'msc': 0.2, 'ace': 0.0, 'als': 0.5, 'eug': 0.3}
sub-141 ses-5 153 0.400 0.350 0.169 0.050 {'msc': 0.0, 'ace': 0.8, 'als': 0.1, 'eug': 0.7}
sub-141 ses-8 198 0.300 0.375 0.174 0.550 {'msc': 0.0, 'ace': 0.2, 'als': 1.0, 'eug': 0.0}
sub-148 ses-2 286 0.250 0.350 0.170 0.075 {'msc': 0.5, 'ace': 0.0, 'als': 0.0, 'eug': 0.5}
sub-148 ses-3 318 0.450 0.375 0.167 0.125 {'msc': 0.9, 'ace': 0.0, 'als': 0.3, 'eug': 0.6}
sub-149 ses-1 219 0.475 0.350 0.172 0.300 {'msc': 0.7, 'ace': 0.0, 'als': 0.3, 'eug': 0.9}
Variant, k at or below the class size (9 other trials of the own odor):
session rois acc_k5 null_mean_k5 null_p99_k5 acc_k9 null_mean_k9 null_p99_k9
sub-134 ses-22 73 0.950 0.222 0.450 0.900 0.203 0.450
sub-141 ses-1 210 0.975 0.220 0.450 0.925 0.206 0.425
sub-141 ses-5 153 1.000 0.219 0.425 1.000 0.207 0.450
sub-141 ses-8 198 0.925 0.217 0.425 0.800 0.206 0.425
sub-148 ses-2 286 0.800 0.213 0.425 0.800 0.203 0.425
sub-148 ses-3 318 0.900 0.215 0.425 0.825 0.202 0.425
sub-149 ses-1 219 0.875 0.209 0.425 0.800 0.204 0.425
Chance 0.25. Main session: the one with 318 ROIs.
Result: partially reproduced. In the 318-glomerulus session (sub-148 ses-3) the 50-s ITI averages classify the preceding odor with 45% accuracy (k = 20), above the 99th percentile of the permutation null (37.5%; null mean 16.7%) but below 50%. Accuracy differs between odors (methyl salicylate 0.9, acetophenone 0.0, allyl sulfide 0.3, eugenol 0.6). The code's choice, cosine distance on the top three PCs, gives 12.5%. Across the seven sessions accuracy is 22.5% to 47.5%, above the null's 99th percentile in three (sub-141 ses-5, sub-148 ses-3, sub-149 ses-1). With k = 20 a trial has 9 neighbours of its own odor among 39, which works against the own odor (null mean 16.7% to 17.4%, below the nominal 25%). With k at or below the class size accuracy is much higher: 90.0% (k = 5) and 82.5% (k = 9) in sub-148 ses-3, and 80.0% to 100% in every session, against null 99th percentiles of 42.5% to 45.0%. All trials of one odor lie in one block, so these accuracies cannot separate the odor from slow differences between blocks (see the pseudo-trials under fig7j).
fig7j-iti-timebin-mouse¶
A time-bin k-NN classifier identifies the preceding odor from glomerular activity throughout the 55 s inter-trial interval; accuracy stays high for its whole duration and falls to chance with shuffled labels.
rng = np.random.default_rng(8)
SEG = 40 # 4-frame bins per 10 s
rows, curves = [], {}
t0 = time.time()
for s in block:
x = olfdyn.dff_session(s["sol"][:, ITI], s) # 40 x 800 x ROIs
n_tr, n_fr, n_roi = x.shape
xb = x.reshape(n_tr, n_fr // 4, 4, n_roi).mean(axis=2) # 40 x 200 x ROIs
X = xb.reshape(-1, n_roi)
trial = np.repeat(np.arange(n_tr), n_fr // 4)
y = np.repeat(s["sol_pres"].odor.to_numpy(), n_fr // 4)
seg = np.tile(np.arange(n_fr // 4) // SEG, n_tr)
nb = olfdyn.knn_neighbours(X, trial, 10, "correlation")
hit = olfdyn.knn_predict(nb, y) == y
null = olfdyn.permuted_accuracy(nb, y, trial, rng, 200, by_bin=seg)
acc_seg = [hit[seg == k].mean() for k in range(5)]
curves[s["name"]] = np.array([hit[np.tile(np.arange(n_fr // 4), n_tr) == k].mean() for k in range(n_fr // 4)])
row = {"session": s["name"], "rois": s["n_rois"]}
for k in range(5):
row[f"acc_{10 * k}-{10 * k + 10}s"] = acc_seg[k]
row[f"p99_{10 * k}-{10 * k + 10}s"] = np.percentile(null[:, k], 99)
row["per_odor_last10s"] = {o: round(float(hit[(seg == 4) & (y == o)].mean()), 2) for o in olfdyn.BLOCK_ODORS}
rows.append(row)
tb = pd.DataFrame(rows)
print(tb.round(3).to_string(index=False))
print(f"(chance 0.25; {time.time() - t0:.0f} s)")
fig, ax = plt.subplots(figsize=(7, 3))
for n, c in curves.items():
ax.plot(np.arange(len(c)) * 0.25, 100 * c, lw=1, label=n)
ax.axhline(25, color="grey", ls=":")
ax.set_xlabel("time after offset (s)"); ax.set_ylabel("accuracy (%)"); ax.set_ylim(0, 100)
ax.legend(fontsize=6, ncol=2)
plt.tight_layout()
print("Plotted accuracies range from", round(100 * min(c.min() for c in curves.values()), 1), "to",
round(100 * max(c.max() for c in curves.values()), 1), "%")
session rois acc_0-10s p99_0-10s acc_10-20s p99_10-20s acc_20-30s p99_20-30s acc_30-40s p99_30-40s acc_40-50s p99_40-50s per_odor_last10s
sub-134 ses-22 73 0.985 0.363 0.991 0.358 0.992 0.342 0.982 0.350 0.973 0.355 {'msc': 1.0, 'ace': 0.94, 'als': 0.96, 'eug': 1.0}
sub-141 ses-1 210 0.984 0.371 0.974 0.378 0.931 0.384 0.919 0.363 0.911 0.388 {'msc': 1.0, 'ace': 0.88, 'als': 0.77, 'eug': 0.99}
sub-141 ses-5 153 0.994 0.351 0.986 0.368 0.976 0.352 0.977 0.356 0.982 0.349 {'msc': 1.0, 'ace': 1.0, 'als': 0.93, 'eug': 1.0}
sub-141 ses-8 198 0.976 0.344 0.954 0.329 0.926 0.330 0.885 0.313 0.884 0.326 {'msc': 1.0, 'ace': 0.98, 'als': 0.72, 'eug': 0.83}
sub-148 ses-2 286 0.932 0.371 0.896 0.350 0.918 0.369 0.935 0.368 0.929 0.364 {'msc': 0.89, 'ace': 0.93, 'als': 0.9, 'eug': 1.0}
sub-148 ses-3 318 0.885 0.335 0.756 0.325 0.755 0.347 0.738 0.359 0.719 0.337 {'msc': 1.0, 'ace': 0.38, 'als': 0.49, 'eug': 1.0}
sub-149 ses-1 219 0.797 0.337 0.732 0.329 0.726 0.333 0.728 0.342 0.722 0.341 {'msc': 0.97, 'ace': 0.3, 'als': 0.66, 'eug': 0.96}
(chance 0.25; 11 s)
Plotted accuracies range from 55.0 to 100.0 %
# Variant: the same classifier on pseudo-trials in the inter-block interval after each solitary block,
# continuing the block's 60.4 s rhythm with no odor given. High accuracy here would show that the
# classifier can name the block from slow differences between blocks, without odor-evoked activity.
rows = []
for s in block:
counts = {o: len(s["pseudo"][o]) for o in olfdyn.BLOCK_ODORS}
use = [o for o in olfdyn.BLOCK_ODORS if counts[o] > 0]
x = np.concatenate([olfdyn.dff_session(s["pseudo"][o], s) for o in use])
n_tr, n_fr, n_roi = x.shape
X = x.reshape(n_tr, n_fr // 4, 4, n_roi).mean(axis=2).reshape(-1, n_roi)
trial = np.repeat(np.arange(n_tr), n_fr // 4)
y = np.repeat(np.concatenate([[o] * counts[o] for o in use]), n_fr // 4)
seg = np.tile(np.arange(n_fr // 4) // SEG, n_tr)
nb = olfdyn.knn_neighbours(X, trial, 10, "correlation")
hit = olfdyn.knn_predict(nb, y) == y
row = {"session": s["name"], "pseudo_trials": counts, "chance": round(max(counts.values()) / n_tr, 3) if use else np.nan}
for k in range(5):
row[f"acc_{10 * k}-{10 * k + 10}s"] = hit[seg == k].mean()
rows.append(row)
print(pd.DataFrame(rows).round(3).to_string(index=False))
print("chance: the share of the largest class")
session pseudo_trials chance acc_0-10s acc_10-20s acc_20-30s acc_30-40s acc_40-50s
sub-134 ses-22 {'msc': 10, 'ace': 10, 'als': 10, 'eug': 10} 0.250 0.944 0.949 0.940 0.958 0.951
sub-141 ses-1 {'msc': 10, 'ace': 10, 'als': 10, 'eug': 10} 0.250 0.899 0.896 0.873 0.895 0.887
sub-141 ses-5 {'msc': 10, 'ace': 10, 'als': 10, 'eug': 10} 0.250 0.938 0.939 0.949 0.953 0.962
sub-141 ses-8 {'msc': 10, 'ace': 10, 'als': 10, 'eug': 10} 0.250 0.919 0.921 0.924 0.911 0.933
sub-148 ses-2 {'msc': 10, 'ace': 10, 'als': 10, 'eug': 7} 0.270 0.955 0.968 0.962 0.963 0.969
sub-148 ses-3 {'msc': 10, 'ace': 10, 'als': 10, 'eug': 9} 0.256 0.698 0.696 0.712 0.703 0.697
sub-149 ses-1 {'msc': 10, 'ace': 10, 'als': 10, 'eug': 8} 0.263 0.746 0.716 0.712 0.715 0.725
chance: the share of the largest class
Result: reproduced. In the 318-glomerulus session the time-bin classifier names the preceding odor with 88.5% accuracy in the first 10 s after offset and 71.9% to 75.6% from 10 to 50 s, each above 50% and above the null's 99th percentile (32.5% to 35.9%). The other six sessions give 72.2% to 99.4% in every segment. This meets the criterion, but the outputs do not show that the accuracy reflects odor memory. The trial-label null does not keep the block structure, and as each odor has one block, no label permutation that keeps whole blocks can serve as a null: renaming whole classes leaves a k-NN classifier's accuracy unchanged except for ties. The same classifier on pseudo-trials in the inter-block intervals, within the 15 min after each block where no odor is given, names the preceding block with 69.6% to 96.9% accuracy, nearly as high as during the ITIs (sub-148 ses-3: 69.6% to 71.2%); in four of seven sessions the pseudo-trial maximum is below the ITI minimum. Slow differences between blocks are enough for most of this accuracy.
figs7e-ibi-chance¶
During the 15 min between blocks odor-specific information is no longer present; classification of 1-min bins of the inter-block interval by the preceding odor is near chance (34.8%).
rng = np.random.default_rng(12)
print("IBI minutes available per session and preceding odor:",
{s["name"]: {o: len(v) for o, v in s["ibi"].items()} for s in block})
rows = []
for s in block:
counts = [len(s["ibi"][o]) for o in olfdyn.BLOCK_ODORS]
X = np.vstack([s["ibi"][o] for o in olfdyn.BLOCK_ODORS]).astype(np.float64)
X = (X - s["F_median"]) / s["F_median"]
y = np.repeat(olfdyn.BLOCK_ODORS, counts)
complete = all(c == 15 for c in counts)
for k in ([50, 20, 10, 5] if complete else [20, 10, 5]):
nb = olfdyn.knn_neighbours(X, np.arange(len(X)), k, "cosine")
pred = olfdyn.knn_predict(nb, y)
rate = np.mean([(pred == y)[y == o].mean() for o in olfdyn.BLOCK_ODORS if (y == o).any()])
null = olfdyn.permuted_accuracy(nb, y, np.arange(len(X)), rng, 1000)
rows.append({"session": s["name"], "complete": complete, "bins": len(X), "k": k, "rate_pct": 100 * rate,
"null_mean_pct": 100 * null.mean(), "null_p95_pct": 100 * np.percentile(null, 95),
"perm_p": float((null >= (pred == y).mean()).mean())})
ib = pd.DataFrame(rows)
print(ib.round(2).to_string(index=False))
for k in [10, 5]:
g = ib[ib.k == k]
print(f"k = {k} (at or below the 14 other bins of the own IBI): rate {g.rate_pct.min():.1f}% to {g.rate_pct.max():.1f}%, "
f"null mean {g.null_mean_pct.min():.1f}% to {g.null_mean_pct.max():.1f}%, permutation p <= {g.perm_p.max():.3f}")
main_ib = ib[ib.complete & (ib.k == 50)]
print(f"\nMain: mean rate over the {len(main_ib)} sessions with four complete IBIs (k=50): {main_ib.rate_pct.mean():.1f}% "
f"(range {main_ib.rate_pct.min():.1f}-{main_ib.rate_pct.max():.1f}); expected 34.8, criterion <= 38.3; chance 25")
IBI minutes available per session and preceding odor: {'sub-134 ses-22': {'msc': 15, 'ace': 15, 'als': 15, 'eug': 15}, 'sub-141 ses-1': {'msc': 15, 'ace': 15, 'als': 15, 'eug': 15}, 'sub-141 ses-5': {'msc': 15, 'ace': 15, 'als': 15, 'eug': 15}, 'sub-141 ses-8': {'msc': 15, 'ace': 15, 'als': 15, 'eug': 15}, 'sub-148 ses-2': {'msc': 15, 'ace': 15, 'als': 15, 'eug': 7}, 'sub-148 ses-3': {'msc': 15, 'ace': 15, 'als': 15, 'eug': 9}, 'sub-149 ses-1': {'msc': 15, 'ace': 15, 'als': 15, 'eug': 8}}
session complete bins k rate_pct null_mean_pct null_p95_pct perm_p sub-134 ses-22 True 60 50 0.00 10.77 23.33 1.00 sub-134 ses-22 True 60 20 90.00 20.62 35.00 0.00 sub-134 ses-22 True 60 10 96.67 22.00 35.00 0.00 sub-134 ses-22 True 60 5 96.67 22.84 35.00 0.00 sub-141 ses-1 True 60 50 13.33 10.91 25.00 0.39 sub-141 ses-1 True 60 20 90.00 19.98 33.33 0.00 sub-141 ses-1 True 60 10 90.00 21.78 35.00 0.00 sub-141 ses-1 True 60 5 98.33 22.70 33.42 0.00 sub-141 ses-5 True 60 50 1.67 10.99 25.00 0.93 sub-141 ses-5 True 60 20 95.00 20.55 33.33 0.00 sub-141 ses-5 True 60 10 98.33 22.10 35.00 0.00 sub-141 ses-5 True 60 5 98.33 22.35 33.33 0.00 sub-141 ses-8 True 60 50 0.00 11.11 26.67 1.00 sub-141 ses-8 True 60 20 73.33 19.86 31.67 0.00 sub-141 ses-8 True 60 10 81.67 22.18 35.00 0.00 sub-141 ses-8 True 60 5 95.00 22.88 35.00 0.00 sub-148 ses-2 False 52 20 66.67 23.66 36.54 0.00 sub-148 ses-2 False 52 10 80.71 25.14 38.56 0.00 sub-148 ses-2 False 52 5 94.76 25.51 38.46 0.00 sub-148 ses-3 False 54 20 35.00 22.34 35.19 0.02 sub-148 ses-3 False 54 10 69.44 23.55 35.19 0.00 sub-148 ses-3 False 54 5 90.00 24.29 35.19 0.00 sub-149 ses-1 False 53 20 58.33 23.53 35.85 0.00 sub-149 ses-1 False 53 10 93.33 24.57 37.74 0.00 sub-149 ses-1 False 53 5 96.67 25.43 39.62 0.00 k = 10 (at or below the 14 other bins of the own IBI): rate 69.4% to 98.3%, null mean 21.8% to 25.1%, permutation p <= 0.000 k = 5 (at or below the 14 other bins of the own IBI): rate 90.0% to 98.3%, null mean 22.4% to 25.5%, permutation p <= 0.000 Main: mean rate over the 4 sessions with four complete IBIs (k=50): 3.8% (range 0.0-13.3); expected 34.8, criterion <= 38.3; chance 25
Result: blocked (methods underspecified). With the Methods' k = 50 on 60 one-minute bins the rate is 0.0% to 13.3% (mean 3.8%), and this procedure's own permutation null averages 10.8%, not the 25% chance of four odors, because the other odors supply most of the 50 neighbours of every bin. The paper's 34.8%, which it calls roughly chance, cannot come from this procedure, so the paper's procedure is not the one its Methods describe. With k at or below the 14 other bins of the own IBI the bins are classified by the preceding block far above chance: 69.4% to 98.3% with k = 10 and 90.0% to 98.3% with k = 5 in all seven sessions (null means 21.8% to 25.5%, permutation p < 0.001), and 35.0% to 95.0% with k = 20. These results contradict "no longer present", but each odor has one IBI per session, so no validation can hold out a block: neighbouring minutes of one IBI are each other's neighbours, and slow drift alone could give these rates.
figs8e-random-order-decline¶
When eleven odors are presented in random order (1 s pulses, 19 s apart), odor information after offset is variable across odors and on average declines to chance after several seconds.
rng = np.random.default_rng(13)
AFTER = slice(48, 48 + 272) # 0 to 17 s after offset (offset 1 s after onset)
rows, curves = [], {}
for s in random_sessions:
q = s["rnd_pres"]
x = olfdyn.dff_session(s["rnd"][:, AFTER], s)
n_tr, n_fr, n_roi = x.shape
xb = x.reshape(n_tr, n_fr // 4, 4, n_roi).mean(axis=2)
X = xb.reshape(-1, n_roi)
trial = np.repeat(np.arange(n_tr), n_fr // 4)
y = np.repeat(q.odor.to_numpy(), n_fr // 4)
tbin = np.tile(np.arange(n_fr // 4), n_tr)
part = np.where(tbin < 8, 0, np.where(tbin >= 48, 2, 1)) # 0-2 s, 2-12 s, 12-17 s
nb = olfdyn.knn_neighbours(X, trial, 10, "correlation")
hit = olfdyn.knn_predict(nb, y) == y
null = olfdyn.permuted_accuracy(nb, y, trial, rng, 200, by_bin=part)
curves[s["name"]] = np.array([hit[tbin == k].mean() for k in range(n_fr // 4)])
early, late = hit[part == 0].mean(), hit[part == 2].mean()
rows.append({"session": s["name"], "rois": s["n_rois"], "rows": len(X), "odors": q.odor.nunique(), "early_0-2s": early,
"null_p99_early": np.percentile(null[:, 0], 99), "late_12-17s": late,
"null_p95_late": np.percentile(null[:, 2], 95),
"early_above": early > np.percentile(null[:, 0], 99), "late_above": late > np.percentile(null[:, 2], 95),
"per_odor_early": {o: round(float(hit[(part == 0) & (y == o)].mean()), 2) for o in sorted(q.odor.unique())}})
rd = pd.DataFrame(rows)
print(rd.drop(columns="per_odor_early").round(3).to_string(index=False))
print("\nPer odor, 0-2 s after offset:")
for _, r in rd.iterrows():
print(" ", r.session, r.per_odor_early)
print(f"\nChance {1 / 11:.3f}. Sessions with early above null: {int(rd.early_above.sum())}; with late not above: {int((~rd.late_above).sum())}; "
f"with both: {int((rd.early_above & ~rd.late_above).sum())} of {len(rd)}")
# The published version's Table 1 (iScience 2026) lists random-order sessions with 175, 340, 323, 306,
# 286, 204, 250 and 280 glomeruli.
table1 = [175, 340, 323, 306, 286, 204, 250, 280]
print("Sessions whose ROI count is in that list:", rd[rd.rois.isin(table1)].session.tolist())
print(rd[rd.rois.isin(table1)][["session", "rois", "early_0-2s", "null_p99_early", "late_12-17s", "null_p95_late", "late_above"]]
.round(3).to_string(index=False))
fig, ax = plt.subplots(figsize=(7, 3))
for n, c in curves.items():
ax.plot(np.arange(len(c)) * 0.25, 100 * c, lw=1, label=n)
ax.axhline(100 / 11, color="grey", ls=":")
ax.set_xlabel("time after offset (s)"); ax.set_ylabel("accuracy (%)"); ax.set_ylim(0, 100); ax.legend(fontsize=7)
plt.tight_layout()
print("Plotted accuracies range from", round(100 * min(c.min() for c in curves.values()), 1), "to",
round(100 * max(c.max() for c in curves.values()), 1), "%")
session rois rows odors early_0-2s null_p99_early late_12-17s null_p95_late early_above late_above
sub-124 ses-3 559 7480 11 0.331 0.107 0.104 0.097 True True
sub-125 ses-1 325 7480 11 0.551 0.123 0.153 0.100 True True
sub-126 ses-1 323 7480 11 0.244 0.106 0.096 0.090 True True
sub-134 ses-6 306 7480 11 0.290 0.116 0.099 0.100 True False
sub-1716 ses-1 403 7480 11 0.568 0.114 0.106 0.094 True True
Per odor, 0-2 s after offset:
sub-124 ses-3 {'2ep': 0.41, 'ace': 0.7, 'als': 0.07, 'ani': 0.15, 'asa': 0.65, 'cin': 0.16, 'eug': 0.19, 'ieg': 0.1, 'mpz': 0.0, 'msc': 0.55, 'tcm': 0.65}
sub-125 ses-1 {'2ep': 0.66, 'ace': 0.69, 'als': 0.29, 'ani': 0.85, 'asa': 0.96, 'cin': 0.14, 'eug': 0.57, 'ieg': 0.05, 'mpz': 0.05, 'msc': 0.8, 'tcm': 1.0}
sub-126 ses-1 {'2ep': 0.26, 'ace': 0.38, 'als': 0.11, 'ani': 0.09, 'asa': 0.29, 'cin': 0.1, 'eug': 0.06, 'ieg': 0.06, 'mpz': 0.1, 'msc': 0.94, 'tcm': 0.3}
sub-134 ses-6 {'2ep': 0.14, 'ace': 0.35, 'als': 0.23, 'ani': 0.42, 'asa': 0.65, 'cin': 0.19, 'eug': 0.33, 'ieg': 0.25, 'mpz': 0.17, 'msc': 0.07, 'tcm': 0.39}
sub-1716 ses-1 {'2ep': 0.45, 'ace': 0.71, 'als': 0.51, 'ani': 0.71, 'asa': 0.9, 'cin': 0.46, 'eug': 0.57, 'ieg': 0.19, 'mpz': 0.1, 'msc': 0.72, 'tcm': 0.91}
Chance 0.091. Sessions with early above null: 5; with late not above: 1; with both: 1 of 5
Sessions whose ROI count is in that list: ['sub-126 ses-1', 'sub-134 ses-6']
session rois early_0-2s null_p99_early late_12-17s null_p95_late late_above
sub-126 ses-1 323 0.244 0.106 0.096 0.09 True
sub-134 ses-6 306 0.290 0.116 0.099 0.10 False
Plotted accuracies range from 4.5 to 83.6 %
Result: partially reproduced. In all five random-order sessions the odor is classified in the first 2 s after offset above the null (24.4% to 56.8% against 99th percentiles of 10.6% to 12.3%; chance 9.1%), with large differences between odors. In the last 5 s accuracy falls to 9.6% to 15.3%, but it is still above the null's 95th percentile in four of the five sessions; only in sub-134 ses-6 is it not. These sessions are mostly not the paper's: the published version's Table 1 lists random-order sessions with 175, 340, 323, 306, 286, 204, 250 and 280 glomeruli, and only sub-126 ses-1 (323) and sub-134 ses-6 (306) have one of these counts. In those two, early accuracy is 24.4% and 29.0%, and late accuracy is above the null's 95th percentile in sub-126 ses-1 (9.6% against 9.0%) and not in sub-134 ses-6 (9.9% against 10.0%). In three sessions odors have different concentration values.
Locust claims (data on Figshare)¶
The locust recordings are not in the dandiset. They are read from the paper's Figshare record, checked against the md5 it gives. Spike counts are in 50 ms bins; trials are on the second axis.
loc2 = olfdyn.figshare_mat("data2_locust_PN66.mat")
loc4 = olfdyn.figshare_mat("data4_locust_PN85.mat")
loc5 = olfdyn.figshare_mat("data5_locust_PN85.mat")
odors66 = list(loc2["odors"])
print("PN66 file: odors", odors66, "array", loc2[odors66[0]].shape, "ON bins", loc2["onBins"][[0, -1]], "OFF bins", loc2["offBins"][[0, -1]])
print("PN85 contrast file: single", list(loc4["singleOdors"]), "pairs", [tuple(r) for r in loc4["contrastOrg"]],
"onInd", loc4["onInd"][[0, -1]], "contrastInd", loc4["contrastInd"][[0, -1]])
print("PN85 persistence file: single", list(loc5["singleOdors"]), "onInd", loc5["onInd"][[0, -1]],
"array", getattr(loc5["data"], "app").shape)
print("Spike counts per bin range:", int(min(loc2[o].min() for o in odors66)), "to", int(max(loc2[o].max() for o in odors66)))
print("PN85 persistence file, variables:", sorted(k for k in loc5 if not k.startswith("__")),
"fields of data:", list(loc5["data"]._fieldnames))
LON, LOFF = slice(200, 280), slice(280, 360) # MATLAB 201:280 and 281:360
PN66 file: odors ['amm', 'bza', 'hex', 'hxa', 'iaa', 'nan', 'noa', 'nen', 'oca', 'phn', 'vny'] array (66, 10, 800) ON bins [201 280] OFF bins [281 360]
PN85 contrast file: single ['app', 'bza', 'cit', 'ger', 'hex', 'iaa', 'oct'] pairs [('app', 'ger', 'appger'), ('app', 'hex', 'apphex'), ('bza', 'ger', 'bzager'), ('bza', 'hex', 'bzahex'), ('cit', 'ger', 'citger'), ('cit', 'hex', 'cithex'), ('iaa', 'ger', 'iaager'), ('iaa', 'hex', 'iaahex'), ('oct', 'ger', 'octger'), ('oct', 'hex', 'octhex')] onInd [301 380] contrastInd [391 470]
PN85 persistence file: single ['app', 'bza', 'cit', 'ger', 'hex', 'iaa', 'oct'] onInd [301 380] array (85, 10, 799)
Spike counts per bin range: 0 to 6
PN85 persistence file, variables: ['contrastInd', 'contrastOdors', 'data', 'nCells', 'nTrials', 'offInd', 'onInd', 'singleOdors'] fields of data: ['app', 'bza', 'cit', 'ger', 'hex', 'iaa', 'oct', 'appger', 'apphex', 'bzager', 'bzahex', 'citger', 'cithex', 'iaager', 'iaahex', 'octger', 'octhex']
fig1g-onoff-angle-locust¶
In the locust antennal lobe the median angle between PN population ON and OFF response vectors exceeds 90 degrees for every one of 11 odors (mean 106.7 +/- 0.84 SEM).
rows = []
for o in odors66:
m = loc2[o].mean(axis=1).T.astype(np.float64) # bins x PNs, trial-averaged
base = m[:200].mean(axis=0)
unc = olfdyn.onoff_angles(m[LON], m[LOFF], center=False)
zero = int((np.linalg.norm(m[LON], axis=1) == 0).sum() + (np.linalg.norm(m[LOFF], axis=1) == 0).sum())
rows.append({"odor": o, **olfdyn.onoff_angles(m[LON], m[LOFF], center=True),
"uncentred": unc["on_off"], "zero_vectors_uncentred": zero,
"ctrl_on_halves": olfdyn.onoff_angles(m[200:240], m[240:280], center=True)["on_off"],
"baseline_subtracted": olfdyn.onoff_angles(m[LON] - base, m[LOFF] - base, center=False)["on_off"]})
la = pd.DataFrame(rows)
print(la.round(2).to_string(index=False))
print(f"\nMain: mean {la.on_off.mean():.2f}, SEM {la.on_off.std(ddof=1) / np.sqrt(len(la)):.2f}, "
f"min {la.on_off.min():.2f}, odors above 90: {int((la.on_off > 90).sum())} of {len(la)}; expected 106.7 +/- 0.84")
print(f"Uncentred: mean {la.uncentred.mean():.2f}; baseline-subtracted: mean {la.baseline_subtracted.mean():.2f}, "
f"above 90: {int((la.baseline_subtracted > 90).sum())}")
print(f"Centring control, first against second half of ON (2 s each): mean {la.ctrl_on_halves.mean():.2f}, "
f"above 90: {int((la.ctrl_on_halves > 90).sum())} of {len(la)}")
odor on_off on_on off_off uncentred zero_vectors_uncentred ctrl_on_halves baseline_subtracted amm 106.52 77.98 70.70 63.14 0 95.97 97.56 bza 104.47 65.65 79.28 65.22 0 97.61 100.33 hex 103.82 78.16 76.77 61.26 0 97.99 92.54 hxa 104.37 73.09 81.15 57.64 0 94.97 93.02 iaa 99.19 79.16 85.85 57.92 0 96.47 93.13 nan 101.08 79.78 78.74 49.09 0 95.03 80.30 noa 108.66 75.23 67.85 54.46 0 95.47 93.45 nen 105.83 73.73 74.79 58.43 0 94.16 91.19 oca 108.11 73.73 71.29 58.03 0 94.06 91.18 phn 107.27 75.51 77.79 64.80 0 99.51 94.98 vny 104.84 75.66 76.66 53.05 0 94.86 83.19 Main: mean 104.92, SEM 0.87, min 99.19, odors above 90: 11 of 11; expected 106.7 +/- 0.84 Uncentred: mean 58.46; baseline-subtracted: mean 91.90, above 90: 9 Centring control, first against second half of ON (2 s each): mean 96.01, above 90: 11 of 11
Result: partially reproduced. All 11 odors give a median ON-OFF angle above 90 degrees (99.2 to 108.7), but the mean, 104.9 degrees (SEM 0.87), is 1.8 degrees below 106.7 +/- 0.84. Without centring, the mean angle is 58.5 degrees; with each PN's pre-odor mean subtracted and no centring it is 91.9, with 9 of 11 above 90. As in the mouse, the centring alone gives angles above 90: the first and second halves of the ON window are 96.0 degrees apart on average (11 of 11 above 90).
figs4a-reference-odor-locust¶
In the locust the ON response to a reference odor is closer to other odors' ON responses than to OFF responses for most reference odors, significantly for 5 of 11.
x = np.vstack([np.vstack([loc2[o].mean(axis=1).T[LON], loc2[o].mean(axis=1).T[LOFF]]) for o in odors66]).astype(np.float64)
xc = x - x.mean(axis=0)
pcs = xc @ np.linalg.svd(xc, full_matrices=False)[2][:3].T
means = np.stack([pcs[i * 80:(i + 1) * 80].mean(0) for i in range(22)])
on, off = means[0::2], means[1::2]
rows = []
for i, o in enumerate(odors66):
a_on = olfdyn.angles(on[i:i + 1], on)[0]
a_off = olfdyn.angles(on[i:i + 1], off)[0]
t = stats.ttest_ind(a_on, a_off)
a_on10 = np.delete(a_on, i) # variant: without the reference's own 0-degree angle
t10 = stats.ttest_ind(a_on10, a_off)
rows.append({"ref": o, "on_mean": a_on.mean(), "off_mean": a_off.mean(), "p": t.pvalue,
"sig_on_smaller": (t.pvalue < 0.05) and a_on.mean() < a_off.mean(),
"on_mean_no_self": a_on10.mean(), "p_no_self": t10.pvalue,
"sig_no_self": (t10.pvalue < 0.05) and a_on10.mean() < a_off.mean()})
rl = pd.DataFrame(rows)
print(rl.round(4).to_string(index=False))
print(f"\nReference odors with ON angles smaller on average: {int((rl.on_mean < rl.off_mean).sum())} of 11; "
f"significant (p<0.05, ON smaller): {int(rl.sig_on_smaller.sum())}; expected 5 +/- 1")
print(f"Variant without the self-angle: ON smaller on average for {int((rl.on_mean_no_self < rl.off_mean).sum())} of 11; "
f"significant: {int(rl.sig_no_self.sum())}")
ref on_mean off_mean p sig_on_smaller on_mean_no_self p_no_self sig_no_self amm 107.3461 52.4769 0.0033 False 118.0807 0.0001 False bza 70.2357 100.5702 0.0520 False 77.2592 0.1016 False hex 64.9514 118.2088 0.0016 True 71.4465 0.0031 True hxa 68.7324 112.3689 0.0085 True 75.6056 0.0163 True iaa 83.0685 79.1280 0.8000 False 91.3754 0.3781 False nan 84.3219 103.2271 0.2535 False 92.7541 0.4757 False noa 87.2613 105.6158 0.2907 False 95.9874 0.5349 False nen 83.7172 91.9828 0.6295 False 92.0889 0.9946 False oca 80.5827 114.7391 0.0592 False 88.6410 0.1155 False phn 73.6846 123.0814 0.0081 True 81.0530 0.0158 True vny 67.4910 124.2046 0.0011 True 74.2401 0.0019 True Reference odors with ON angles smaller on average: 9 of 11; significant (p<0.05, ON smaller): 4; expected 5 +/- 1 Variant without the self-angle: ON smaller on average for 8 of 11; significant: 4
Result: reproduced. Four reference odors have ON angles significantly smaller than OFF angles (hexanol, hexanoic acid, phenethylamine, 4-vinyl anisole), within 5 +/- 1, and 9 of 11 have smaller ON angles on average. Ammonia differs significantly in the other direction (ON 107.3, OFF 52.5 degrees, p = 0.0033). The 11 ON angles include the reference's own angle of 0 degrees, as the criterion states; without it the count is still 4 and ON is smaller on average for 8 of 11. Benzaldehyde (p = 0.052) and octanoic acid (p = 0.059) are near the threshold.
fig5e-contrast-locust¶
In the locust, PN responses to geraniol presented 0.5 s after another odor are less correlated with that odor's solitary response than solitary geraniol is (one-tailed t-test p < 0.05, 5 odor pairs).
d4 = loc4["data"]
on_i, ct_i = loc4["onInd"] - 1, loc4["contrastInd"] - 1
def vec(name, idx):
return getattr(d4, name).astype(np.float64).mean(axis=1)[:, idx].mean(axis=1)
for second in ["ger", "hex"]:
rows = []
for first in ["app", "bza", "cit", "iaa", "oct"]:
rows.append({"first": first, "r_solitary": olfdyn.corr(vec(first, on_i), vec(second, on_i)),
"r_sequential": olfdyn.corr(vec(first, on_i), vec(first + second, ct_i))})
lc = pd.DataFrame(rows)
t = stats.ttest_rel(lc.r_solitary, lc.r_sequential, alternative="greater")
t2 = stats.ttest_rel(lc.r_solitary, lc.r_sequential)
print(f"--- second odor {second}{' (main)' if second == 'ger' else ' (variant)'}")
print(lc.round(3).to_string(index=False))
print(f"solitary {lc.r_solitary.mean():.3f}, sequential {lc.r_sequential.mean():.3f}, one-tailed paired p = {t.pvalue:.3g} "
f"(two-tailed {t2.pvalue:.3g}, as MATLAB's ttest default in Script_Fig4.m)")
--- second odor ger (main) first r_solitary r_sequential app 0.584 0.485 bza 0.196 0.055 cit 0.638 0.294 iaa 0.447 0.127 oct 0.565 0.177 solitary 0.486, sequential 0.227, one-tailed paired p = 0.00553 (two-tailed 0.0111, as MATLAB's ttest default in Script_Fig4.m) --- second odor hex (variant) first r_solitary r_sequential app 0.406 0.406 bza 0.106 0.213 cit 0.354 0.252 iaa 0.575 0.309 oct 0.592 0.427 solitary 0.407, sequential 0.321, one-tailed paired p = 0.128 (two-tailed 0.257, as MATLAB's ttest default in Script_Fig4.m)
Result: reproduced. Solitary first odors correlate with solitary geraniol at r = 0.486 on average and with geraniol presented after them at 0.227, lower for all five first odors (one-tailed paired p = 0.0055). With hexanol as the second odor the difference is smaller (0.407 against 0.321, p = 0.13).
fig7f-iti-knn-locust¶
In the locust, PN activity averaged over the 20 s after odor offset classifies the preceding odor (k-NN on angular distance), so ITI activity is odor specific.
rng = np.random.default_rng(17)
d5 = loc5["data"]
single = list(loc5["singleOdors"])
X = np.vstack([getattr(d5, o).astype(np.float64)[:, :, 380:780].mean(axis=2).T for o in single]) # 70 x 85
y = np.repeat(single, 10)
zero = int((np.linalg.norm(X, axis=1) == 0).sum())
nb = olfdyn.knn_neighbours(X, np.arange(len(X)), 10, "cosine")
acc = (olfdyn.knn_predict(nb, y) == y).mean()
null = olfdyn.permuted_accuracy(nb, y, np.arange(len(X)), rng, 1000)
print(f"Main (k=10, full-dimensional angular distance, 20 s): accuracy {acc:.3f}, null mean {null.mean():.3f}, "
f"99th percentile {np.percentile(null, 99):.3f}; chance {1 / 7:.3f}; zero vectors {zero}")
print("Per odor:", {o: round(float((olfdyn.knn_predict(nb, y) == y)[y == o].mean()), 2) for o in single})
X21 = np.vstack([getattr(d5, o).astype(np.float64)[:, :, 380:799].mean(axis=2).T for o in single])
xc = X21 - X21.mean(0)
pcs = xc @ np.linalg.svd(xc, full_matrices=False)[2][:3].T
nb3 = olfdyn.knn_neighbours(pcs, np.arange(len(X)), 4, "cosine")
print(f"Variant (code: k=4, cosine on top three PCs, 21 s): accuracy {(olfdyn.knn_predict(nb3, y) == y).mean():.3f}")
Main (k=10, full-dimensional angular distance, 20 s): accuracy 1.000, null mean 0.114, 99th percentile 0.257; chance 0.143; zero vectors 0
Per odor: {'app': 1.0, 'bza': 1.0, 'cit': 1.0, 'ger': 1.0, 'hex': 1.0, 'iaa': 1.0, 'oct': 1.0}
Variant (code: k=4, cosine on top three PCs, 21 s): accuracy 0.971
Result: reproduced. The 20-s ITI averages classify all 70 trials correctly (100%; null 99th percentile 25.7%, chance 14.3%). The code's choice gives 97.1%. The published version's Table 1 gives this dataset a solitary block structure, so all trials of one odor were likely given together, as in the mouse. The Figshare file holds only the 799 bins (40 s) around each trial, so no pseudo-trial control between blocks can be built from it, and these outputs cannot separate the odor from slow differences between blocks.
fig7i-iti-timebin-locust¶
In the locust a time-bin k-NN classifier's accuracy decays across the 20 s after offset to about 40% averaged over odors, but for most odors does not fall to chance.
rng = np.random.default_rng(19)
rows_X, y, trial, tbin = [], [], [], []
for i, o in enumerate(single):
a = getattr(d5, o).astype(np.float64)[:, :, 380:799] # PNs x trials x 419 bins
nb_ = int(np.ceil(a.shape[2] / 4))
b = np.stack([a[:, :, j * 4:(j + 1) * 4].mean(axis=2) for j in range(nb_)], axis=2) # last bin partial, as the code
for tr in range(10):
rows_X.append(b[:, tr, :].T)
y += [o] * nb_; trial += [i * 10 + tr] * nb_; tbin += list(range(nb_))
X = np.vstack(rows_X); y = np.array(y); trial = np.array(trial); tbin = np.array(tbin)
const = X.std(axis=1) == 0
print(f"Time bins: {len(X)}; bins with no spikes or constant counts (dropped from queries and neighbours): {int(const.sum())}")
Xv, yv, trv, tbv = X[~const], y[~const], trial[~const], tbin[~const]
nb = olfdyn.knn_neighbours(Xv, trv, 10, "correlation")
hit = olfdyn.knn_predict(nb, yv) == yv
last = tbv >= tbv.max() - 24 # last 25 bins, 5 s
per_odor = {o: hit[last & (yv == o)].mean() for o in single}
null_rows = []
trials_u = np.unique(trv)
lab = np.array([yv[trv == t][0] for t in trials_u])
pos = np.searchsorted(trials_u, trv)
for _ in range(200):
perm = rng.permutation(lab)[pos]
h = olfdyn.knn_predict(nb, perm) == perm
null_rows.append([h[last & (perm == o)].mean() for o in single] + [h[last].mean()])
null_rows = np.array(null_rows)
print(f"Mean over odors in the last 5 s: {np.mean(list(per_odor.values())):.3f} (expected 0.40 +/- 0.04); chance {1 / 7:.3f}")
print(f"Over the first 2 s after offset: {hit[tbv < 10].mean():.3f}; over 0-21 s: {hit.mean():.3f}")
print("Per odor in the last 5 s, with the 99th percentile of the null:")
for j, o in enumerate(single):
print(f" {o}: {per_odor[o]:.3f} (null p99 {np.percentile(null_rows[:, j], 99):.3f})")
print("Odors above the null:", sum(per_odor[o] > np.percentile(null_rows[:, j], 99) for j, o in enumerate(single)), "of 7")
curve = np.array([hit[tbv == k].mean() for k in range(tbv.max() + 1)])
fig, ax = plt.subplots(figsize=(6, 2.5))
ax.plot(np.arange(len(curve)) * 0.2, 100 * curve, color="k")
ax.axhline(100 / 7, color="grey", ls=":"); ax.set_ylim(0, 100)
ax.set_xlabel("time after offset (s)"); ax.set_ylabel("accuracy (%)")
plt.tight_layout()
print("Plotted accuracies range from", round(100 * curve.min(), 1), "to", round(100 * curve.max(), 1), "%")
Time bins: 7350; bins with no spikes or constant counts (dropped from queries and neighbours): 0
Mean over odors in the last 5 s: 0.290 (expected 0.40 +/- 0.04); chance 0.143 Over the first 2 s after offset: 0.787; over 0-21 s: 0.396 Per odor in the last 5 s, with the 99th percentile of the null: app: 0.380 (null p99 0.276) bza: 0.244 (null p99 0.228) cit: 0.228 (null p99 0.200) ger: 0.264 (null p99 0.176) hex: 0.540 (null p99 0.168) iaa: 0.160 (null p99 0.148) oct: 0.216 (null p99 0.140) Odors above the null: 7 of 7 Plotted accuracies range from 18.6 to 94.3 %
Result: partially reproduced. Accuracy decays from 78.7% in the first 2 s after offset to 29.0% averaged over odors in the last 5 s, below 40 +/- 4. In the last 5 s every one of the 7 odors stays above its null's 99th percentile (16.0% to 54.0% against 14.0% to 27.6%), so for most odors accuracy does not fall to chance. Averaged over the whole 0 to 21 s the accuracy is 39.6%, a reading of "about 40%" the criterion did not list. As for fig7f, the files allow no control for the block structure.
Summary¶
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.
| Claim | Data | Outcome |
|---|---|---|
| mouse-sessions | dandiset | reproduced |
| figshare-trace-link | with inputs | reproduced |
| fig2h-onoff-angle-mouse | dandiset | partially reproduced |
| fig2f-onoff-anticorrelated-mouse | dandiset | not reproduced |
| fig4b-reference-odor-mouse | dandiset | partially reproduced |
| fig4d-clustering-mouse | dandiset | not reproduced |
| fig6e-contrast-msc | dandiset | not reproduced |
| fig6e-contrast-ace | dandiset | blocked: data missing |
| fig7d-persistent-anticorrelation | dandiset | not reproduced |
| fig7h-iti-knn-mouse | dandiset | partially reproduced |
| fig7j-iti-timebin-mouse | dandiset | reproduced |
| figs7e-ibi-chance | dandiset | blocked: methods underspecified |
| figs8e-random-order-decline | dandiset | partially reproduced |
| figs2-pid-decay | dandiset | blocked: data missing |
| fig1g-onoff-angle-locust | elsewhere | partially reproduced |
| figs4a-reference-odor-locust | elsewhere | reproduced |
| fig5e-contrast-locust | elsewhere | reproduced |
| fig7f-iti-knn-locust | elsewhere | reproduced |
| fig7i-iti-timebin-locust | elsewhere | partially reproduced |
How the analysis differs from the paper's.
- dF/F is computed from the dandiset's raw F with F0 from the 2 s before each onset (responses) or the session median (ITI and IBI). The paper's deposited arrays are close to F z-scored within each 57-s trial window. With that z-scoring the contrast enhancement (fig6e) and the negative correlation of later activity with the ON response from 5 s after offset (fig7d) appear, and neither appears with dF/F or with F z-scored on the 2 s before onset. The ON/OFF anticorrelation of fig2f does not appear with any of these (median r 0.068 with the whole-window z-score, 12 of 28 below 0). On synthetic glomeruli with excitation only, the whole-window z-score gives negative correlations between ON and later activity by itself. The z-scored results are shown as variants and do not decide the outcomes.
- The ON-OFF angles follow the code's centring, which makes the mean ON and OFF vectors opposite (uncentred, the mouse mean is 71.8 degrees; in the locust, with the pre-odor baseline subtracted and no centring, it is 91.9, with 9 of 11 odors above 90); the same procedure gives about 101 degrees between the two halves of the mouse ON window and 96 between those of the locust ON window, so angles above 90 here do not show inversion.
- Windows are placed in seconds; the paper's code counts frames at a nominal 16 per second, which lasts 51.2 s to 66.4 s for a nominal 57-s trial (912 frames) at the sessions' 13.7 to 17.8 Hz, so the windows drift by 6 to 9 s at the extremes (about 1.8 s in the 15.5-Hz session of the deposited array). Frame-counted windows are shown as a variant for fig2h, fig2f and fig7d and give the same conclusions.
- Sessions are matched to figures by ROI count. The contrast claim uses all seven block sessions, as the preprint does not identify its three; the three that the published version's Table 1 names are shown as a variant. The random-order claim uses five sessions, of which two have ROI counts that appear in that table.
- The k-NN choices follow the Methods where they can be run; the Methods' k = 50 cannot be run on the 40 ITI averages and is degenerate on the 60 IBI bins. With k at or below the class size, the ITI averages and the IBI bins are both classified at 69% to 100%. In the block sessions every trial of one odor lies in one block, and pseudo-trials within the 15 min after a block, with no odor, are classified nearly as well as the ITIs, so these classifiers cannot separate the odor from slow differences between blocks.
What was found about the dandiset.
- The seven block-protocol sessions in four mice match the paper's mouse dataset, and the raw F of the 198-glomerulus session matches the paper's deposited array (median within-trial r 0.968, same ROI order).
- The fluorescence is raw F (unit "n.a.") with no description of its extraction or of the normalisation the paper used, which decides several outcomes here.
- The trials table holds 1-s valve rows without a description; blocks, solitary and sequential presentations must be rebuilt from timing. Some sessions fill the gaps with "Empty" rows and some do not, and concentrations have no unit and are 0.0 for two odorants in three sessions.
- In three files the raw two-photon movie holds half the frames of the traces.
- No sequential pair ends in acetophenone, and there are no PID recordings, so Figure 6E (right), S6 and S2 cannot be tested. The random-order session of Figure S8D/E is not in the dandiset; five other random-order sessions are, four of them from mice not in the block sessions, with odor names that differ from the paper's for two odors.
- Genotype and age are not stated, and every session starts at 00:00.
- The dandiset puts the 318-ROI session under sub-148 with the 286-ROI session; the published version's Table 1 puts it under mouse 4 with the 219-ROI session (sub-149 here).
print(f"Notebook ran in {(time.time() - started) / 60:.1f} min")
Notebook ran in 2.2 min