Talbot et al. (2018) on CRCNS hc-16: single cells, coordination and frame-specific firing¶
Paper: Talbot, Sparks, Dvorak, Curran, Alarcon and Fenton (2018), "Normal CA1 Place Fields but Discoordinated Network Discharge in a Fmr1-Null Mouse Model of Fragile X Syndrome", Neuron 97:684-697, doi:10.1016/j.neuron.2017.12.043.
Dataset: CRCNS hc-16 (Sparks, Talbot, Dvorak and Fenton), doi:10.6080/K01V5C44, read from the
public bucket crcnsarchive at hc-16/live/data.zip, streamed with fsspec and zipfile.
Version: CRCNS data sets are not versioned, so the version is recorded as unversioned. The
paper cites the data set's DOI only (Key Resources Table), and the bucket holds one copy of the
data, under live/ (there is no old/), so that copy is the one read. The cell below reads the
bucket's checksums.md5; the file read is data.zip, MD5 45b534dd6dbeb94db539bc429223eecb
(data descriptor version 0.6, April 26, 2018). A later change of the zip would change that sum.
Claims: claims/crcns-hc-16/10.1016-j.neuron.2017.12.043.yaml, written from the version of
record and its supplement (Figures S1 to S6, Tables S1 and S2), as its opening comment says.
| Claim | Here |
|---|---|
| fig1ab-baseline-transmission, fig1cd-potentiation, figs1-avoidance-equal | blocked: slice recordings and slice-cohort behaviour not in the data set |
| unit-and-cell-counts | blocked: the data set holds 15 of the paper's 21 mice and no isolation measures (counts printed for the record) |
| interneuron-proportion | tested (rate and bursting only; no waveforms) |
| place-cell-proportion | tested |
| figs3-nonspatial-properties-similar | tested (no waveform width) |
| fig2-spatial-quality-similar | tested (large box only) |
| tables1-place-cell-properties-similar | tested |
| tables1-nonspatial-active-pixels | tested |
| figs4-remapping-similar | blocked: no box and cylinder recordings of the same cells (shown below) |
| fig3b-place-cell-modulation-similar, fig3bc-weaker-modulation-nonspatial-interneurons, fig3b-phase-preference, fig3d-stereotypy, figs5a-power-spectra-similar | blocked: no LFP in the data set |
| overdispersion-passes, overdispersion-cells | tested |
| fig4b-covariance-30-passes, fig4b-covariance-40-passes | tested |
| fig4cd-kendall-place-cells, fig4d-network-state-recurrence, kendall-nonspatial-pyramidal | tested |
| kendall-quiet-states | tested (low-speed periods; no sleep recordings) |
| locomotion-similar | tested |
| fig5b-frame-specific-ipos, fig5c-ensemble-coordination, figs6-room-frame-preference | tested |
| fig6-covariance-sem-predicts-errors, fig6-kendall-predicts-errors, discussion-single-cell-measures-no-relation, fig6-avoidance-behavior | blocked: avoidance sessions not in the data set |
Population used throughout: the data set has no isolation measures (IsoIBG, IsoINN) and no waveforms, so the paper's well-isolated cells cannot be selected. Every cluster other than cluster 0 of a tetrode is used. This is the largest difference from the paper's method and bears on every tested claim. Two variants bear on it: a refractory filter (clusters with at most 0.5% of interspike intervals below 2 ms), a proxy for isolation, and single-cell measures on the code's 32 x 32 grid, which needs no spatial scale.
In [1]:
import sys, time, hashlib, re
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy import stats
import hc16
pd.set_option("display.width", 200)
pd.set_option("display.max_columns", 40)
pd.set_option("display.max_rows", 200)
T0 = time.time()
md5 = hc16.check_md5_list()
print("checksums.md5:", md5)
assert md5["data.zip"] == hc16.DATA_ZIP_MD5
z = hc16.open_zip()
files = hc16.list_members(z)
sessions = hc16.load_sessions(files)
S = {s.id: s for s in sessions}
print(f"loaded {len(sessions)} sessions in {time.time() - T0:.0f} s")
checksums.md5: {'data.zip': '45b534dd6dbeb94db539bc429223eecb', 'docs/crcns_hc-16_data_description.docx': '41bc3a8b453215e070ad5b925b3fe27b', 'docs/crcns_hc-16_data_description.pdf': '86544fa66ee16226829909e6ec996715', 'filelist.txt': '8d60384041e6f525496910b5e6b52fa2'}
loaded 41 sessions in 29 s
In [2]:
print("Zip members matched by name:", len(files), "of", sum(1 for i in z.infolist() if not i.is_dir()),
"(the rest are .DS_Store)")
print(files.groupby(["env", "genotype", "frame", "part", "times"]).size().rename("files").to_string())
sess = pd.DataFrame([dict(session=s.id, env=s.env, genotype=s.genotype, mouse=s.mouse, date=s.date,
cond=s.cond, cut=s.cut, rate_hz=round(s.rate, 3), duration_s=round(s.duration, 1),
frames="+".join(sorted(s.frames)), clusters=len(s.counts), cells=len(s.cells),
untracked=round(1 - s.tracked.mean(), 3),
extent_px=f"{s.extent_px[0]:.0f}x{s.extent_px[1]:.0f}",
px_per_cm=round(s.px_per_cm, 2), leftover_bytes=s.leftover)
for s in sessions])
print(sess.to_string())
print("\nOne full TS header (%s, first file):" % sessions[0].id)
print(sessions[0].header)
print("Header fields over every file read:", sorted(set().union(*[s.header_fields for s in sessions])))
print("SCALE_X and SCALE_Y values in the first file of each session:",
sorted({v for s in sessions for v in re.findall(r"%SCALE_[XY]\S*\s+([\d.]+)", s.header)}))
print("\nMice per genotype:", sess.groupby("genotype").mouse.nunique().to_dict())
print("Sessions per environment and genotype:\n", sess.groupby(["env", "genotype"]).size().unstack())
print("Position sampling rate (Hz): %.3f to %.3f; durations %.0f to %.0f s"
% (sess.rate_hz.min(), sess.rate_hz.max(), sess.duration_s.min(), sess.duration_s.max()))
# cluster 0 of each tetrode
c0 = [(s.id, u, s.counts[u].sum() / s.duration) for s in sessions for u in s.counts if u.endswith("C0")]
c0r = np.array([r for _, _, r in c0])
print(f"\nCluster 0: {len(c0)} files, rate median {np.median(c0r):.1f} AP/s (range {c0r.min():.1f} to {c0r.max():.1f})")
allr = np.array([s.counts[u].sum() / s.duration for s in sessions for u in s.cells])
print(f"Other clusters: {len(allr)}, rate median {np.median(allr):.2f} AP/s (range {allr.min():.3f} to {allr.max():.1f})")
# spikes are sorted, inside the session, and each lies near its position sample
bad_sort = 0; offs = []
for s in sessions:
for u in s.cells:
sp = s.spikes[u]
bad_sort += int(np.any(np.diff(sp) < 0))
k = np.repeat(np.arange(s.nsamp), s.counts[u])
offs.append(((sp - s.t0) / 1e4 * s.rate - k))
offs = np.concatenate(offs)
print(f"Clusters with unsorted spike times: {bad_sort}")
print("Spike time minus the nominal time of its position sample (t0 + k/rate), in samples: "
"percentiles 0/50/100 = %.2f / %.2f / %.2f" % tuple(np.percentile(offs, [0, 50, 100])))
res = []
for s in sessions:
for u in s.counts:
k = np.repeat(np.arange(s.nsamp), s.counts[u])
res.append((s.spikes[u] / 1e4 - s.sample_time[k]) * s.rate)
res = np.concatenate(res)
print("Spike time minus the estimated time of its sample (hc16.Session.sample_time), in samples: "
"percentiles 0.1/50/99.9 = %.2f / %.2f / %.2f" % tuple(np.percentile(res, [0.1, 50, 99.9])))
for s in sessions:
for u in s.cells:
d = np.diff(s.spikes[u])
if np.any(d < 0):
print(f" unsorted: {s.id} {u}: {int((d < 0).sum())} backward steps of {len(d) + 1} spikes, "
f"largest {-d.min() / 10:.1f} ms")
# x and y range (tracker pixels, 0 to 255; (0, 0) is untracked)
for env, g in sess.groupby("env"):
print(f"{env}: room-frame extent {', '.join(g.extent_px)} px; px/cm {g.px_per_cm.min()} to {g.px_per_cm.max()}")
Zip members matched by name: 4087 of 4090 (the rest are .DS_Store)
env genotype frame part times
big_box Fmr1-null Rm full False 175
True 175
pos full False 7
True 7
WT Rm _0-15min False 4
full False 139
True 139
pos _0-15min False 2
full False 8
True 8
rotating_cylinder Fmr1-null Ar _0-15min False 65
_15-30min False 40
full False 97
True 97
Rm _0-15min False 65
_15-30min False 40
full False 97
True 97
pos _0-15min False 3
_15-30min False 2
full False 4
True 4
WT Ar _0-15min False 94
_15-30min False 94
full False 94
True 94
Rm _0-15min False 94
_15-30min False 94
full False 94
True 94
pos _0-15min False 3
_15-30min False 3
full False 3
True 3
stationary_cylinder Fmr1-null Ar _0-15min False 43
_15-30min False 43
full False 43
True 43
Rm _0-15min False 221
_15-30min False 221
full False 285
True 285
pos _0-15min False 8
_15-30min False 8
full False 10
True 10
WT Rm _0-15min False 141
_15-30min False 141
full False 257
True 257
pos _0-15min False 7
_15-30min False 7
full False 9
True 9
session env genotype mouse date cond cut rate_hz duration_s frames clusters cells untracked extent_px px_per_cm leftover_bytes
0 M10_052214_blackbox big_box Fmr1-null M10 052214 blackbox cl120915 30.019 1798.6 Rm 6 5 0.020 216x220 4.36 0
1 M18_100814_blackbox big_box Fmr1-null M18 100814 blackbox cl120915 30.019 1798.4 Rm 19 16 0.060 186x182 3.68 0
2 M18_100914_rotating rotating_cylinder Fmr1-null M18 100914 rotating cl120815 30.009 1799.2 Ar+Rm 32 24 0.045 231x234 5.81 0
3 M18_100914_stationary1 stationary_cylinder Fmr1-null M18 100914 stationary1 cl120815 30.009 1799.1 Rm 32 24 0.049 232x234 5.83 0
4 M18_100914_stationary2 stationary_cylinder Fmr1-null M18 100914 stationary2 cl120815 30.009 1799.2 Rm 32 24 0.046 228x232 5.75 0
5 M29_061215_blackbox1 big_box Fmr1-null M29 061215 blackbox1 cl121115 30.019 1199.3 Rm 32 24 0.031 190x181 3.71 0
6 M29_112414_rotating rotating_cylinder Fmr1-null M29 112414 rotating cl020316 30.019 1798.0 Ar+Rm 25 17 0.038 220x225 5.56 0
7 M29_112414_stationary1 stationary_cylinder Fmr1-null M29 112414 stationary1 cl120715 30.019 1798.7 Ar+Rm 25 17 0.048 227x228 5.68 0
8 M29_112414_stationary2 stationary_cylinder Fmr1-null M29 112414 stationary2 cl020316 30.018 1798.7 Rm 25 17 0.027 226x227 5.67 0
9 M33_051715_blackbox big_box Fmr1-null M33 051715 blackbox cl121515 30.019 1798.8 Rm 25 18 0.050 186x181 3.67 0
10 M33_051815_rotating rotating_cylinder Fmr1-null M33 051815 rotating cl112715 30.009 1799.5 Ar+Rm 9 2 0.013 234x226 5.75 0
11 M33_051815_stationary1 stationary_cylinder Fmr1-null M33 051815 stationary1 cl112715 30.010 1799.4 Ar+Rm 9 2 0.023 234x233 5.84 0
12 M33_051815_stationary2 stationary_cylinder Fmr1-null M33 051815 stationary2 cl112715 30.009 1781.8 Ar+Rm 9 2 0.068 228x225 5.66 0
13 M34_041015_blackbox big_box Fmr1-null M34 041015 blackbox cl121615 30.019 1798.4 Rm 42 34 0.057 188x182 3.70 0
14 M34_041115_stationary1 stationary_cylinder Fmr1-null M34 041115 stationary1 cl101115 30.009 1798.7 Rm 46 38 0.044 235x235 5.88 0
15 M34_041115_stationary2 stationary_cylinder Fmr1-null M34 041115 stationary2 cl101115b 30.009 1798.9 Rm 46 38 0.040 234x234 5.85 0
16 M37_050815_blackbox big_box Fmr1-null M37 050815 blackbox cl121615 30.019 1798.9 Rm 34 26 0.059 186x180 3.66 0
17 M37_051115_rotating rotating_cylinder Fmr1-null M37 051115 rotating cl101015 30.009 1799.5 Ar+Rm 31 23 0.029 226x226 5.65 0
18 M37_051115_stationary1 stationary_cylinder Fmr1-null M37 051115 stationary1 cl101015 30.009 1799.5 Rm 30 22 0.030 228x231 5.74 0
19 M37_051115_stationary2 stationary_cylinder Fmr1-null M37 051115 stationary2 cl101015 30.009 1799.4 Rm 31 23 0.038 227x230 5.71 0
20 M9_052614_blackbox big_box Fmr1-null M9 052614 blackbox cl121515 30.020 1798.2 Rm 17 12 0.021 219x221 4.40 0
21 M1_040714_blackbox big_box WT M1 040714 blackbox cl120915 30.020 1798.3 Rm 10 9 0.004 188x207 3.95 0
22 M11_071414_blackbox big_box WT M11 071414 blackbox cl123115A 30.019 1798.6 Rm 10 8 0.032 218x222 4.40 0
23 M16_090414_rotating rotating_cylinder WT M16 090414 rotating cl101015 30.009 1799.1 Ar+Rm 58 50 0.084 235x235 5.88 0
24 M16_090414_stationary1 stationary_cylinder WT M16 090414 stationary1 cl101015 30.009 1799.1 Rm 58 50 0.051 231x233 5.80 0
25 M16_090414_stationary2 stationary_cylinder WT M16 090414 stationary2 cl101015 30.009 1799.0 Rm 58 50 0.046 232x232 5.80 0
26 M16_091214_blackbox big_box WT M16 091214 blackbox cl123115 30.019 1198.9 Rm 35 27 0.013 184x180 3.64 0
27 M28_032515_rotating rotating_cylinder WT M28 032515 rotating cl020216 30.011 1792.2 Ar+Rm 13 10 0.116 231x230 5.76 0
28 M28_032515_stationary1 stationary_cylinder WT M28 032515 stationary1 cl020216 30.009 1799.0 Rm 13 10 0.129 222x225 5.59 0
29 M28_032515_stationary2 stationary_cylinder WT M28 032515 stationary2 cl020216 30.009 1798.8 Rm 18 10 0.139 227x229 5.69 0
30 M28_060315_blackbox1 big_box WT M28 060315 blackbox1 cl123115 30.019 1198.4 Rm 15 14 0.050 187x179 3.66 0
31 M3_041314_blackbox big_box WT M3 041314 blackbox cl123115 30.019 1798.3 Rm 5 3 0.026 177x192 3.69 0
32 M32_041715_rotating rotating_cylinder WT M32 041715 rotating cl101015 30.008 1799.0 Ar+Rm 23 15 0.049 223x229 5.65 0
33 M32_041715_stationary1 stationary_cylinder WT M32 041715 stationary1 cl101015 30.009 1799.2 Rm 23 15 0.053 223x222 5.56 0
34 M32_041715_stationary2 stationary_cylinder WT M32 041715 stationary2 cl101015 30.009 1799.1 Rm 23 15 0.043 214x223 5.47 0
35 M32_060315_blackbox1 big_box WT M32 060315 blackbox1 cl123115 30.018 1198.6 Rm 32 24 0.048 182x178 3.60 0
36 M35_040115_blackbox big_box WT M35 040115 blackbox cl123115 30.019 1798.7 Rm 23 15 0.012 189x180 3.69 0
37 M35_040215_stationary1 stationary_cylinder WT M35 040215 stationary1 cl100915 30.009 1799.1 Rm 21 13 0.028 227x230 5.71 0
38 M35_040215_stationary2 stationary_cylinder WT M35 040215 stationary2 cl100915 30.009 1798.4 Rm 21 13 0.020 224x226 5.62 0
39 M36_041415_blackbox big_box WT M36 041415 blackbox cl123115 30.019 1798.6 Rm 9 1 0.022 185x181 3.66 0
40 M36_050115_stationary1 stationary_cylinder WT M36 050115 stationary1 cl120215 30.009 1799.0 Rm 22 14 0.026 226x228 5.67 0
One full TS header (M10_052214_blackbox, first file):
13
%%BEGIN_HEADER
%SAMPLING_INTERVAL(samps/sec) 30.019213
%FIRST_TIMESTAMP(100usUnits) 7430
%SCALE_Y(RatioTracktoMapPixels) 1.00
%SCALE_X(RatioTracktoMapPixels) 1.00
TIME_SERIES FORMAT:
Ascii header:
NumberOfLinesInHeader
(NumberOfLinesInHeader - 1) of Ascii Info.
Binary data format:
Xposition:1byte Yposition:1byte NspksInSample:1byte TimestampOfSpk-1:4bytes...TimeStampOfSpk-n:4bytes
%%ENDHEADER
Header fields over every file read: ['BEGIN_HEADER', 'ENDHEADER', 'FIRST_TIMESTAMP', 'SAMPLING_INTERVAL', 'SCALE_X', 'SCALE_Y']
SCALE_X and SCALE_Y values in the first file of each session: ['1.00']
Mice per genotype: {'Fmr1-null': 7, 'WT': 8}
Sessions per environment and genotype:
genotype Fmr1-null WT
env
big_box 7 8
rotating_cylinder 4 3
stationary_cylinder 10 9
Position sampling rate (Hz): 30.008 to 30.020; durations 1198 to 1800 s
Cluster 0: 273 files, rate median 25.9 AP/s (range 0.0 to 184.7)
Other clusters: 774, rate median 1.39 AP/s (range 0.001 to 52.2)
Clusters with unsorted spike times: 6 Spike time minus the nominal time of its position sample (t0 + k/rate), in samples: percentiles 0/50/100 = -0.85 / 9.54 / 118.20
Spike time minus the estimated time of its sample (hc16.Session.sample_time), in samples: percentiles 0.1/50/99.9 = -0.01 / 0.50 / 1.01 unsorted: M34_041015_blackbox P0C2: 1 backward steps of 21439 spikes, largest 830.1 ms unsorted: M34_041015_blackbox P2C2: 1 backward steps of 6603 spikes, largest 217.9 ms unsorted: M37_051115_stationary1 P2C2: 1 backward steps of 50269 spikes, largest 1441.8 ms unsorted: M37_051115_stationary1 P8C2: 1 backward steps of 88956 spikes, largest 1129.8 ms unsorted: M16_090414_rotating P1C1: 1 backward steps of 35444 spikes, largest 604.9 ms unsorted: M35_040215_stationary2 P0C2: 1 backward steps of 47283 spikes, largest 399.9 ms big_box: room-frame extent 216x220, 186x182, 190x181, 186x181, 188x182, 186x180, 219x221, 188x207, 218x222, 184x180, 187x179, 177x192, 182x178, 189x180, 185x181 px; px/cm 3.6 to 4.4 rotating_cylinder: room-frame extent 231x234, 220x225, 234x226, 226x226, 235x235, 231x230, 223x229 px; px/cm 5.56 to 5.88 stationary_cylinder: room-frame extent 232x234, 228x232, 227x228, 226x227, 234x233, 228x225, 235x235, 234x234, 228x231, 227x230, 231x233, 232x232, 222x225, 227x229, 223x222, 214x223, 227x230, 224x226, 226x228 px; px/cm 5.47 to 5.88
Checks of the data against the description, and of the files against each other:
- The descriptor gives x and y in 0-255 pixels and times in 100 us units, but no spatial scale. The scale is taken from the extent of the tracked room-frame positions (0.1 to 99.9 percentiles) and the enclosure sizes of the STAR Methods (large box 50 cm, arena 40 cm). This is an assumption; the printed px/cm differ between box sessions (about 3.6 and 4.4), which fits two camera set-ups and is consistent with, not proof of, the assumed sizes.
- The TS headers hold the sampling rate, the first timestamp and a track-to-map scale of 1.00; no isolation measure and no spatial scale in cm (printed above).
- The descriptor's example header has a sampling interval of 58.3 samples/s; every file here has 30.0, matching the 30 frames/s of the STAR Methods.
- The header's sampling rate and first timestamp do not give the times of the position samples: spikes lie up to 118 samples (about 4 s) after the nominal time t0 + k/rate of the sample they are filed under, as the rate drifts against the spike clock and frames are dropped. The helper estimates each sample's time from the spikes filed under it; after that the spikes lie 0 to 1 sample after their sample's time, as they should. Six spike trains have one backward step in their times (printed above).
- The next cells check that the
.timesfiles repeat the spike times of the TS files, that the_0-15minand_15-30minfiles are parts of the full sessions (one position sample at the split is in neither part, so a spike filed under it is lost), that the arena-frame files of stationary sessions repeat the room frame, and whether any recording is deposited twice.
In [3]:
rng = np.random.default_rng(0)
# .times files against the TS spikes, for every cell of three sessions
chk = []
for sid in ["M16_090414_stationary1", "M34_041015_blackbox", "M37_051115_rotating"]:
s = S[sid]
for u in s.cells:
nm = files[(files.mouse == s.mouse) & (files.date == s.date) & (files.cond == s.cond) & (files.unit == u)
& (files.frame == "Rm") & (files.part == "full") & files.times].name.iloc[0]
v = np.array(z.read(nm).split(), dtype=np.int64)
chk.append(np.array_equal(v, s.spikes[u]))
print(f".times files identical to the TS spike times: {sum(chk)} of {len(chk)} cells (3 sessions)")
# 0-15 and 15-30 min parts against the full session
s = S["M16_090414_rotating"]
part = hc16.read_members([f"data/rotating_cylinder/WT/M16_090414_rotating-cl101015-fix_0-15min.{u}Rm" for u in s.cells[:5]]
+ [f"data/rotating_cylinder/WT/M16_090414_rotating-cl101015-fix_15-30min.{u}Rm" for u in s.cells[:5]])
for u in s.cells[:5]:
a = part[f"data/rotating_cylinder/WT/M16_090414_rotating-cl101015-fix_0-15min.{u}Rm"]
b = part[f"data/rotating_cylinder/WT/M16_090414_rotating-cl101015-fix_15-30min.{u}Rm"]
joined = np.r_[a["spk"], b["spk"]]
print(f"{u}: 0-15 min {len(a['x'])} samples, 15-30 min {len(b['x'])}, full {s.nsamp}; "
f"spikes {len(a['spk'])}+{len(b['spk'])} vs {len(s.spikes[u])}, "
f"parts equal to the full spike train: {np.array_equal(joined, s.spikes[u])}")
na, nb_ = len(a["x"]), len(b["x"])
print("First part equals the full file's first samples:", np.array_equal(a["x"], s.frames["Rm"][0][:na]),
"; second part equals its last samples:", np.array_equal(b["x"], s.frames["Rm"][0][s.nsamp - nb_:]))
dropped = np.arange(na, s.nsamp - nb_)
print("Samples of the full file in neither part:", dropped.tolist(),
"; spikes there per cell:", {u: int(s.counts[u][dropped].sum()) for u in s.cells[:5]})
# arena frame of stationary sessions
for s in sessions:
if "Ar" in s.frames and s.env == "stationary_cylinder":
print(f"{s.id}: arena-frame positions identical to room frame: "
f"{np.array_equal(s.frames['Ar'][0], s.frames['Rm'][0]) and np.array_equal(s.frames['Ar'][1], s.frames['Rm'][1])}")
for s in sessions:
if s.env == "rotating_cylinder":
same = np.mean((s.frames["Ar"][0] == s.frames["Rm"][0]) & (s.frames["Ar"][1] == s.frames["Rm"][1]))
print(f"{s.id}: fraction of samples with equal room and arena positions {same:.3f}")
# recordings deposited twice: identical spike trains or identical position series across sessions
h_spk, h_pos = {}, {}
for s in sessions:
h_pos.setdefault(hashlib.md5(s.frames["Rm"][0].tobytes() + s.frames["Rm"][1].tobytes()).hexdigest(), []).append(s.id)
for u in s.cells:
if len(s.spikes[u]) > 10:
h_spk.setdefault(hashlib.md5(s.spikes[u].tobytes()).hexdigest(), []).append(f"{s.id}.{u}")
dup_pos = [v for v in h_pos.values() if len(v) > 1]
dup_spk = [v for v in h_spk.values() if len(v) > 1]
dup_spk_x = [v for v in dup_spk if len({x.split('.')[0] for x in v}) > 1]
print(f"Position series shared by more than one session: {dup_pos}")
print(f"Identical spike trains across sessions: {len(dup_spk_x)}; within a session: {len(dup_spk) - len(dup_spk_x)}")
print(dup_spk[:10])
.times files identical to the TS spike times: 107 of 107 cells (3 sessions)
P0C1: 0-15 min 26981 samples, 15-30 min 27007, full 53989; spikes 20173+17268 vs 37441, parts equal to the full spike train: True
P0C2: 0-15 min 26981 samples, 15-30 min 27007, full 53989; spikes 12203+11069 vs 23273, parts equal to the full spike train: False
P0C3: 0-15 min 26981 samples, 15-30 min 27007, full 53989; spikes 587+558 vs 1145, parts equal to the full spike train: True
P0C4: 0-15 min 26981 samples, 15-30 min 27007, full 53989; spikes 650+483 vs 1133, parts equal to the full spike train: True
P0C5: 0-15 min 26981 samples, 15-30 min 27007, full 53989; spikes 743+634 vs 1377, parts equal to the full spike train: True
First part equals the full file's first samples: True ; second part equals its last samples: True
Samples of the full file in neither part: [26981] ; spikes there per cell: {'P0C1': 0, 'P0C2': 1, 'P0C3': 0, 'P0C4': 0, 'P0C5': 0}
M29_112414_stationary1: arena-frame positions identical to room frame: True
M33_051815_stationary1: arena-frame positions identical to room frame: True
M33_051815_stationary2: arena-frame positions identical to room frame: True
M18_100914_rotating: fraction of samples with equal room and arena positions 0.048
M29_112414_rotating: fraction of samples with equal room and arena positions 0.041
M33_051815_rotating: fraction of samples with equal room and arena positions 0.015
M37_051115_rotating: fraction of samples with equal room and arena positions 0.030
M16_090414_rotating: fraction of samples with equal room and arena positions 0.067
M28_032515_rotating: fraction of samples with equal room and arena positions 0.118
M32_041715_rotating: fraction of samples with equal room and arena positions 0.051
Position series shared by more than one session: []
Identical spike trains across sessions: 0; within a session: 0
[]
Clock check. Spikes and positions are in the same file, each spike filed under a position sample. To show that spike counts and positions share a clock, the spatial information of every cluster is computed with the true alignment and with the spike counts shifted in time against the positions; if the alignment were off, the true alignment would not stand out.
In [4]:
def info_shifted(s, shift_s):
idx, shape = hc16.pixel_index(s, "Rm")
occ = hc16.occupancy(s, idx, shape)
k = int(round(shift_s * s.rate))
out = []
for u in s.cells:
n = s.counts[u]
s.counts[u] = np.roll(n, k)
rm, o = hc16.rate_map(s, u, idx, shape, occ)
s.counts[u] = n
out.append(hc16.information(rm, o))
return np.array(out)
shifts = [0, 0.5, 2, 10, 60, 300]
rows = []
for s in sessions:
if s.env == "rotating_cylinder":
continue
for sh in shifts:
v = info_shifted(s, sh)
rows.append(dict(session=s.id, shift_s=sh, median_info=np.nanmedian(v)))
clk = pd.DataFrame(rows).pivot(index="session", columns="shift_s", values="median_info")
print("Median information content (bits/spike) of each session's clusters, by shift of spikes against positions:")
print(clk.round(3).to_string())
print("Sessions where the true alignment gives the highest median:", int((clk.idxmax(axis=1) == 0).sum()), "of", len(clk))
Median information content (bits/spike) of each session's clusters, by shift of spikes against positions: shift_s 0.0 0.5 2.0 10.0 60.0 300.0 session M10_052214_blackbox 1.402 1.394 1.287 0.954 0.895 0.899 M11_071414_blackbox 1.106 0.874 0.568 0.464 0.474 0.470 M16_090414_stationary1 0.845 0.858 0.691 0.583 0.598 0.561 M16_090414_stationary2 0.894 0.853 0.754 0.663 0.650 0.666 M16_091214_blackbox 2.607 2.395 2.177 2.168 2.188 2.236 M18_100814_blackbox 2.616 2.613 2.607 2.470 2.381 2.477 M18_100914_stationary1 0.673 0.668 0.609 0.547 0.552 0.540 M18_100914_stationary2 0.537 0.515 0.471 0.462 0.450 0.474 M1_040714_blackbox 1.762 1.734 1.523 1.214 1.267 1.224 M28_032515_stationary1 2.747 2.751 2.715 2.692 2.696 2.675 M28_032515_stationary2 1.676 1.633 1.484 1.236 1.283 1.214 M28_060315_blackbox1 1.307 1.304 1.145 0.947 0.994 0.990 M29_061215_blackbox1 1.120 1.084 1.017 0.989 0.968 0.954 M29_112414_stationary1 0.165 0.176 0.153 0.126 0.122 0.124 M29_112414_stationary2 0.302 0.275 0.190 0.178 0.164 0.167 M32_041715_stationary1 1.950 1.929 1.234 0.996 0.999 0.963 M32_041715_stationary2 2.002 1.973 1.496 1.307 1.323 1.423 M32_060315_blackbox1 1.731 1.769 1.480 1.283 1.230 1.268 M33_051715_blackbox 3.172 3.105 2.862 2.720 2.924 2.866 M33_051815_stationary1 1.578 1.556 1.197 1.083 1.083 1.077 M33_051815_stationary2 1.993 1.808 1.934 1.907 2.179 2.357 M34_041015_blackbox 1.109 1.102 1.012 0.934 0.899 0.911 M34_041115_stationary1 0.314 0.294 0.259 0.216 0.206 0.211 M34_041115_stationary2 0.311 0.301 0.231 0.181 0.178 0.173 M35_040115_blackbox 2.158 2.087 1.987 1.772 1.694 1.616 M35_040215_stationary1 0.455 0.422 0.329 0.317 0.340 0.313 M35_040215_stationary2 0.448 0.427 0.282 0.292 0.291 0.281 M36_041415_blackbox 1.176 1.055 0.790 0.683 0.707 0.723 M36_050115_stationary1 0.927 0.917 0.758 0.683 0.679 0.670 M37_050815_blackbox 2.094 1.939 1.838 1.579 1.678 1.700 M37_051115_stationary1 0.454 0.417 0.347 0.319 0.307 0.297 M37_051115_stationary2 0.440 0.381 0.329 0.296 0.293 0.297 M3_041314_blackbox 1.825 1.737 1.516 1.593 1.519 1.800 M9_052614_blackbox 0.738 0.695 0.570 0.460 0.472 0.474 Sessions where the true alignment gives the highest median: 29 of 34
In 29 of 34 stationary sessions the true alignment gives the highest median information content, and the median falls with the shift in most of the rest; this is consistent with spikes and positions sharing a clock, not proof of it.
Cell table and classification¶
One row per cluster and session, with the single-cell measures on unsmoothed 1.25 cm room-frame maps. Classes (STAR Methods, "Classification of Cell Types"; the waveform criterion cannot be applied): interneuron, rate above 2 AP/s; pyramidal, rate below 2 AP/s with an ISI histogram peak below 10 ms; place cell, a pyramidal cell with rate 0.1-2 AP/s, coherence above 0.4 and information above 0.4 bits/spike ("Characterization of Place Cells"). Clusters below 2 AP/s without bursts are left unclassified. A variant drops the burst criterion (rate alone).
"Standard recordings" (the session a cell is classified in) is open in the paper; the claims
file names two readings. R1: each cell once, from its first stationary-environment session (the
large box, or the first stationary arena session of a day). R2: every stationary session (large
box, stationary arena 1 and 2) on its own, a cell counted per session. Sessions of a day share a
cluster cutting file (the cl... part of the name), so a cluster name within a cut is taken to be
the same cell across that day's sessions (an assumption: the files carry no other identifier).
In [5]:
cells = hc16.cell_table(sessions)
cells["cls"] = hc16.classify(cells)
cells["cls_rate"] = hc16.classify(cells, use_burst=False)
stationary = cells[cells.env != "rotating_cylinder"].copy()
first_sid = {}
for s in sessions:
if s.env == "rotating_cylinder":
continue
first_sid.setdefault(f"{s.mouse}_{s.date}", s.id) # sessions are sorted: stationary1 before stationary2
R1 = stationary[stationary.session.isin(set(first_sid.values()))].copy()
R2 = stationary.copy()
print("R1 sessions:", sorted(first_sid.values()))
for name, t in [("R1", R1), ("R2", R2)]:
print(f"\n{name}: classes by genotype (rate and bursting)")
print(t.groupby(["genotype", "cls"]).size().unstack(fill_value=0))
print(f"{name}: classes by genotype (rate only)")
print(t.groupby(["genotype", "cls_rate"]).size().unstack(fill_value=0))
print("\nRanges of the measures used (all stationary sessions):")
print(R2[["rate", "coherence", "info", "active", "burst_prob", "isi_mean", "isi_peak_ms", "infield_mean", "infield_20pk"]]
.describe().round(3).to_string())
R1 sessions: ['M10_052214_blackbox', 'M11_071414_blackbox', 'M16_090414_stationary1', 'M16_091214_blackbox', 'M18_100814_blackbox', 'M18_100914_stationary1', 'M1_040714_blackbox', 'M28_032515_stationary1', 'M28_060315_blackbox1', 'M29_061215_blackbox1', 'M29_112414_stationary1', 'M32_041715_stationary1', 'M32_060315_blackbox1', 'M33_051715_blackbox', 'M33_051815_stationary1', 'M34_041015_blackbox', 'M34_041115_stationary1', 'M35_040115_blackbox', 'M35_040215_stationary1', 'M36_041415_blackbox', 'M36_050115_stationary1', 'M37_050815_blackbox', 'M37_051115_stationary1', 'M3_041314_blackbox', 'M9_052614_blackbox']
R1: classes by genotype (rate and bursting)
cls interneuron nonspatial place unclassified
genotype
Fmr1-null 106 48 51 33
WT 68 60 48 27
R1: classes by genotype (rate only)
cls_rate interneuron nonspatial place
genotype
Fmr1-null 106 78 54
WT 68 84 51
R2: classes by genotype (rate and bursting)
cls interneuron nonspatial place unclassified
genotype
Fmr1-null 166 62 70 44
WT 103 91 55 42
R2: classes by genotype (rate only)
cls_rate interneuron nonspatial place
genotype
Fmr1-null 166 103 73
WT 103 129 59
Ranges of the measures used (all stationary sessions):
rate coherence info active burst_prob isi_mean isi_peak_ms infield_mean infield_20pk
count 633.000 633.000 633.000 633.000 632.000 632.000 619.000 558.000 280.000
mean 4.622 0.377 1.560 0.528 0.142 10.103 13.006 8.113 14.850
std 7.771 0.283 1.645 0.346 0.112 37.935 28.797 9.591 9.660
min 0.001 -0.067 0.011 0.001 0.000 0.020 1.500 0.254 1.099
25% 0.332 0.154 0.256 0.204 0.045 0.215 4.500 2.129 8.306
50% 1.382 0.338 1.024 0.536 0.131 0.709 5.500 4.601 13.037
75% 4.581 0.536 2.250 0.890 0.215 2.971 8.500 9.503 19.166
max 50.305 1.335 8.968 0.997 0.519 509.042 393.500 56.872 53.765
Two variants of the population and of the maps, shown with the claims they bear on. The data set gives no isolation measure, so the fraction of interspike intervals below 2 ms (refractory violations) is shown as a proxy: a cluster mixing several units has many. The refractory filter keeps clusters with at most 0.5% of intervals below 2 ms, a standard choice. The spatial scale is an assumption, so the single-cell measures are also computed on the code's 32 x 32 grid of the 256 x 256 tracker frame (8 tracker pixels a side), which needs no scale.
In [6]:
print("Fraction of ISIs below 2 ms by class, all sessions:")
print(cells.groupby("cls").refractory.describe(percentiles=[0.5, 0.9]).round(4).to_string())
ISO = 0.005
print(f"Clusters with more than {ISO:.1%} of ISIs below 2 ms:", int((cells.refractory > ISO).sum()), "of", len(cells))
cells_code = hc16.cell_table(sessions, index_fn=hc16.pixel_index_code)
cells_code["cls"] = hc16.classify(cells_code)
print("32 x 32 grid pixel side in cm at the assumed scale:",
{env: f"{8 / g.px_per_cm.max():.2f} to {8 / g.px_per_cm.min():.2f}" for env, g in sess.groupby("env")})
print("\nClass on 1.25 cm pixels (rows) against the 32 x 32 grid (columns), all sessions:")
print(pd.crosstab(cells.cls, cells_code.cls).to_string())
R1_code = cells_code.loc[R1.index]; R2_code = cells_code.loc[R2.index]
R1_iso = R1[R1.refractory <= ISO]; R2_iso = R2[R2.refractory <= ISO]
for name, t in [("R1, 32 x 32 grid", R1_code), ("R2, 32 x 32 grid", R2_code),
("R1, refractory filter", R1_iso), ("R2, refractory filter", R2_iso)]:
print(f"\n{name}: classes by genotype")
print(t.groupby(["genotype", "cls"]).size().unstack(fill_value=0).to_string())
print("\nInterneuron rate (AP/s), R1: quantiles by genotype")
print(R1[R1.cls == "interneuron"].groupby("genotype").rate.describe(percentiles=[0.25, 0.5, 0.75]).round(2).to_string())
Fraction of ISIs below 2 ms by class, all sessions:
count mean std min 50% 90% max
cls
interneuron 328.0 0.0028 0.0049 0.0 0.0008 0.0071 0.0339
nonspatial 189.0 0.0013 0.0066 0.0 0.0000 0.0024 0.0769
place 141.0 0.0006 0.0014 0.0 0.0000 0.0018 0.0104
unclassified 116.0 0.0005 0.0023 0.0 0.0000 0.0007 0.0211
Clusters with more than 0.5% of ISIs below 2 ms: 70 of 774
32 x 32 grid pixel side in cm at the assumed scale: {'big_box': '1.82 to 2.22', 'rotating_cylinder': '1.36 to 1.44', 'stationary_cylinder': '1.36 to 1.46'}
Class on 1.25 cm pixels (rows) against the 32 x 32 grid (columns), all sessions:
cls interneuron nonspatial place unclassified
cls
interneuron 328 0 0 0
nonspatial 0 156 33 0
place 0 12 129 0
unclassified 0 0 0 116
R1, 32 x 32 grid: classes by genotype
cls interneuron nonspatial place unclassified
genotype
Fmr1-null 106 36 63 33
WT 68 53 55 27
R2, 32 x 32 grid: classes by genotype
cls interneuron nonspatial place unclassified
genotype
Fmr1-null 166 51 81 44
WT 103 81 65 42
R1, refractory filter: classes by genotype
cls interneuron nonspatial place unclassified
genotype
Fmr1-null 91 47 51 33
WT 64 57 47 26
R2, refractory filter: classes by genotype
cls interneuron nonspatial place unclassified
genotype
Fmr1-null 134 59 70 44
WT 95 86 54 41
Interneuron rate (AP/s), R1: quantiles by genotype
count mean std min 25% 50% 75% max
genotype
Fmr1-null 106.0 11.47 10.67 2.08 3.54 6.95 16.60 50.30
WT 68.0 7.27 7.04 2.03 2.83 4.34 8.47 34.48
unit-and-cell-counts¶
Of 1,115 clusters, 499 came from 12 WT and 616 from 9 Fmr1-null mice; 226 WT and 334 Fmr1-null well-isolated cells formed the data set.
Criterion: mice counted exactly (11 or 12 WT); clusters and cells counted once or per session. The data set holds no isolation measures, so the well-isolated cells cannot be counted; the cluster counts are given with and without cluster 0.
In [7]:
mice = sess.groupby("genotype").mouse.unique()
print("Mice in the data set:", {g: sorted(v) for g, v in mice.items()})
print("Mice per genotype:", sess.groupby("genotype").mouse.nunique().to_dict(), "(paper: WT 12, Fmr1-null 9)")
cnt = []
for s in sessions:
cnt.append(dict(session=s.id, genotype=s.genotype, env=s.env, mouse=s.mouse, day=f"{s.mouse}_{s.date}",
cut=s.cut, clusters=len(s.cells), with_c0=len(s.counts)))
cnt = pd.DataFrame(cnt)
# once: each day's cut counted once (a cluster name within a cut is one cell across the day's sessions)
once = cnt.drop_duplicates(["day", "cut"])
once_day = cnt.sort_values("session").drop_duplicates(["day"])
tab = pd.DataFrame({
"once per cut file": once.groupby("genotype").clusters.sum(),
"once per day (R1)": R1.groupby("genotype").size(),
"per session, stationary (R2)": R2.groupby("genotype").size(),
"per session, all sessions": cnt.groupby("genotype").clusters.sum(),
"once per cut, with cluster 0": once.groupby("genotype").with_c0.sum(),
"per session, all, with cluster 0": cnt.groupby("genotype").with_c0.sum(),
})
print("\nClusters (cluster 0 left out unless said):")
print(tab.T.to_string())
print("Paper: clusters WT 499, Fmr1-null 616; well-isolated cells WT 226, Fmr1-null 334 (tolerance 10%)")
for col in tab.columns:
w, k = tab.loc["WT", col], tab.loc["Fmr1-null", col]
print(f" {col}: WT {w} ({'in' if abs(w - 499) <= 49.9 else 'out of'} 499 +/- 49.9), "
f"Fmr1-null {k} ({'in' if abs(k - 616) <= 61.6 else 'out of'} 616 +/- 61.6); "
f"as cells: WT {'in' if abs(w - 226) <= 22.6 else 'out of'} 226 +/- 22.6, "
f"Fmr1-null {'in' if abs(k - 334) <= 33.4 else 'out of'} 334 +/- 33.4")
Mice in the data set: {'Fmr1-null': ['M10', 'M18', 'M29', 'M33', 'M34', 'M37', 'M9'], 'WT': ['M1', 'M11', 'M16', 'M28', 'M3', 'M32', 'M35', 'M36']}
Mice per genotype: {'Fmr1-null': 7, 'WT': 8} (paper: WT 12, Fmr1-null 9)
Clusters (cluster 0 left out unless said):
genotype Fmr1-null WT
once per cut file 294 203
once per day (R1) 238 203
per session, stationary (R2) 342 291
per session, all sessions 408 366
once per cut, with cluster 0 389 276
per session, all, with cluster 0 557 490
Paper: clusters WT 499, Fmr1-null 616; well-isolated cells WT 226, Fmr1-null 334 (tolerance 10%)
once per cut file: WT 203 (out of 499 +/- 49.9), Fmr1-null 294 (out of 616 +/- 61.6); as cells: WT out of 226 +/- 22.6, Fmr1-null out of 334 +/- 33.4
once per day (R1): WT 203 (out of 499 +/- 49.9), Fmr1-null 238 (out of 616 +/- 61.6); as cells: WT out of 226 +/- 22.6, Fmr1-null out of 334 +/- 33.4
per session, stationary (R2): WT 291 (out of 499 +/- 49.9), Fmr1-null 342 (out of 616 +/- 61.6); as cells: WT out of 226 +/- 22.6, Fmr1-null in 334 +/- 33.4
per session, all sessions: WT 366 (out of 499 +/- 49.9), Fmr1-null 408 (out of 616 +/- 61.6); as cells: WT out of 226 +/- 22.6, Fmr1-null out of 334 +/- 33.4
once per cut, with cluster 0: WT 276 (out of 499 +/- 49.9), Fmr1-null 389 (out of 616 +/- 61.6); as cells: WT out of 226 +/- 22.6, Fmr1-null out of 334 +/- 33.4
per session, all, with cluster 0: WT 490 (in 499 +/- 49.9), Fmr1-null 557 (in 616 +/- 61.6); as cells: WT out of 226 +/- 22.6, Fmr1-null out of 334 +/- 33.4
Result: blocked (data missing). The data set holds 8 WT and 7 Fmr1-null mice, where the paper has 12 (or 11) and 9: it is a subset of the paper's data, so its counts cannot test the paper's. The well-isolated cells cannot be selected either, as the data set has no isolation measures. For the record, the cluster counts are within 10% of 499 and 616 only when every session is counted on its own and cluster 0 (the unsorted spikes of each tetrode) is included, 490 WT and 557 Fmr1-null, which with 15 of 21 mice is a coincidence; counted once per cut file without cluster 0 there are 203 WT and 294 Fmr1-null clusters.
figs4-remapping-similar (blocked)¶
The claim compares the same cells in a square box and the circular arena. The table pairs each mouse's box sessions with its arena sessions, with their dates and cut files.
In [8]:
rows = []
for m, g in sess.groupby("mouse"):
for _, b in g[g.env == "big_box"].iterrows():
for _, a in g[g.env != "big_box"].iterrows():
rows.append(dict(mouse=m, genotype=b.genotype, box=b.session, box_cut=b.cut, arena=a.session,
arena_cut=a.cut, same_day=b.date == a.date, same_cut=b.cut == a.cut))
fx = pd.DataFrame(rows)
print(fx.to_string())
print("Box and arena session pairs:", len(fx), "; on the same day:", int(fx.same_day.sum()),
"; with the same cut file:", int(fx.same_cut.sum()))
print("Small box sessions in the data set:", int((files.env == "small_box").sum()))
mouse genotype box box_cut arena arena_cut same_day same_cut 0 M16 WT M16_091214_blackbox cl123115 M16_090414_rotating cl101015 False False 1 M16 WT M16_091214_blackbox cl123115 M16_090414_stationary1 cl101015 False False 2 M16 WT M16_091214_blackbox cl123115 M16_090414_stationary2 cl101015 False False 3 M18 Fmr1-null M18_100814_blackbox cl120915 M18_100914_rotating cl120815 False False 4 M18 Fmr1-null M18_100814_blackbox cl120915 M18_100914_stationary1 cl120815 False False 5 M18 Fmr1-null M18_100814_blackbox cl120915 M18_100914_stationary2 cl120815 False False 6 M28 WT M28_060315_blackbox1 cl123115 M28_032515_rotating cl020216 False False 7 M28 WT M28_060315_blackbox1 cl123115 M28_032515_stationary1 cl020216 False False 8 M28 WT M28_060315_blackbox1 cl123115 M28_032515_stationary2 cl020216 False False 9 M29 Fmr1-null M29_061215_blackbox1 cl121115 M29_112414_rotating cl020316 False False 10 M29 Fmr1-null M29_061215_blackbox1 cl121115 M29_112414_stationary1 cl120715 False False 11 M29 Fmr1-null M29_061215_blackbox1 cl121115 M29_112414_stationary2 cl020316 False False 12 M32 WT M32_060315_blackbox1 cl123115 M32_041715_rotating cl101015 False False 13 M32 WT M32_060315_blackbox1 cl123115 M32_041715_stationary1 cl101015 False False 14 M32 WT M32_060315_blackbox1 cl123115 M32_041715_stationary2 cl101015 False False 15 M33 Fmr1-null M33_051715_blackbox cl121515 M33_051815_rotating cl112715 False False 16 M33 Fmr1-null M33_051715_blackbox cl121515 M33_051815_stationary1 cl112715 False False 17 M33 Fmr1-null M33_051715_blackbox cl121515 M33_051815_stationary2 cl112715 False False 18 M34 Fmr1-null M34_041015_blackbox cl121615 M34_041115_stationary1 cl101115 False False 19 M34 Fmr1-null M34_041015_blackbox cl121615 M34_041115_stationary2 cl101115b False False 20 M35 WT M35_040115_blackbox cl123115 M35_040215_stationary1 cl100915 False False 21 M35 WT M35_040115_blackbox cl123115 M35_040215_stationary2 cl100915 False False 22 M36 WT M36_041415_blackbox cl123115 M36_050115_stationary1 cl120215 False False 23 M37 Fmr1-null M37_050815_blackbox cl121615 M37_051115_rotating cl101015 False False 24 M37 Fmr1-null M37_050815_blackbox cl121615 M37_051115_stationary1 cl101015 False False 25 M37 Fmr1-null M37_050815_blackbox cl121615 M37_051115_stationary2 cl101015 False False Box and arena session pairs: 26 ; on the same day: 0 ; with the same cut file: 0 Small box sessions in the data set: 0
No box session shares a day or a cut file with an arena session of the same mouse, and the data set has no waveforms, so no cell can be followed from the square to the circular enclosure. The claim is blocked (data missing): the small box sessions recorded on the days of arena sessions (Table S2, for example M32 on 04-17-15) are not in the data set.
interneuron-proportion¶
A lower proportion of the well-isolated cells is classified as interneurons in WT (28, 12.4%) than in Fmr1-null mice (98, 29.3%).
Criterion: percentages within tolerance (12.4 +/- 1.24, 29.3 +/- 2.93) and a two-proportion z test with p < 0.05, WT lower. Classification by rate and bursting, the variant the claims file names when waveforms are missing; rate alone is shown too. Readings R1 and R2 for the session.
In [9]:
def proportions(t, col, cls, denom_classes=None):
out = {}
for g in ["WT", "Fmr1-null"]:
tt = t[t.genotype == g]
d = tt if denom_classes is None else tt[tt[col].isin(denom_classes)]
out[g] = (int((tt[col] == cls).sum()), len(d))
z, p = hc16.prop_test(out["WT"][0], out["WT"][1], out["Fmr1-null"][0], out["Fmr1-null"][1])
return out, z, p
rows = []
for rname, t in [("R1", R1), ("R2", R2), ("R1 refractory filter", R1_iso), ("R2 refractory filter", R2_iso)]:
for col in ["cls", "cls_rate"]:
for den_name, den in [("all clusters", None), ("classified cells", ["interneuron", "place", "nonspatial"])]:
if col == "cls_rate" and den is not None:
continue
o, zz, pp = proportions(t, col, "interneuron", den)
rows.append(dict(reading=rname, classes="rate+burst" if col == "cls" else "rate only",
denominator=den_name, WT=f"{o['WT'][0]}/{o['WT'][1]}",
WT_pct=100 * o["WT"][0] / o["WT"][1], KO=f"{o['Fmr1-null'][0]}/{o['Fmr1-null'][1]}",
KO_pct=100 * o["Fmr1-null"][0] / o["Fmr1-null"][1], z=zz, p_two_sided=pp))
ip = pd.DataFrame(rows)
ip["WT_in_tol"] = (ip.WT_pct - 12.4).abs() <= 1.24
ip["KO_in_tol"] = (ip.KO_pct - 29.3).abs() <= 2.93
print(ip.round(4).to_string())
print("\nBy mouse (R1, rate+burst, all clusters): interneuron fraction")
print(R1.groupby(["genotype", "mouse"]).apply(lambda d: pd.Series(dict(n=len(d), interneurons=(d.cls == "interneuron").sum(),
frac=(d.cls == "interneuron").mean()))).round(3).to_string())
reading classes denominator WT WT_pct KO KO_pct z p_two_sided WT_in_tol KO_in_tol
0 R1 rate+burst all clusters 68/203 33.4975 106/238 44.5378 -2.3643 0.0181 False False
1 R1 rate+burst classified cells 68/176 38.6364 106/205 51.7073 -2.5535 0.0107 False False
2 R1 rate only all clusters 68/203 33.4975 106/238 44.5378 -2.3643 0.0181 False False
3 R2 rate+burst all clusters 103/291 35.3952 166/342 48.5380 -3.3337 0.0009 False False
4 R2 rate+burst classified cells 103/249 41.3655 166/298 55.7047 -3.3406 0.0008 False False
5 R2 rate only all clusters 103/291 35.3952 166/342 48.5380 -3.3337 0.0009 False False
6 R1 refractory filter rate+burst all clusters 64/194 32.9897 91/222 40.9910 -1.6838 0.0922 False False
7 R1 refractory filter rate+burst classified cells 64/168 38.0952 91/189 48.1481 -1.9128 0.0558 False False
8 R1 refractory filter rate only all clusters 64/194 32.9897 91/222 40.9910 -1.6838 0.0922 False False
9 R2 refractory filter rate+burst all clusters 95/276 34.4203 134/307 43.6482 -2.2779 0.0227 False False
10 R2 refractory filter rate+burst classified cells 95/235 40.4255 134/263 50.9506 -2.3526 0.0186 False False
11 R2 refractory filter rate only all clusters 95/276 34.4203 134/307 43.6482 -2.2779 0.0227 False False
By mouse (R1, rate+burst, all clusters): interneuron fraction
n interneurons frac
genotype mouse
Fmr1-null M10 5.0 2.0 0.400
M18 40.0 13.0 0.325
M29 41.0 22.0 0.537
M33 20.0 4.0 0.200
M34 72.0 40.0 0.556
M37 48.0 18.0 0.375
M9 12.0 7.0 0.583
WT M1 9.0 4.0 0.444
M11 8.0 5.0 0.625
M16 77.0 24.0 0.312
M28 24.0 7.0 0.292
M3 3.0 1.0 0.333
M32 39.0 12.0 0.308
M35 28.0 10.0 0.357
M36 15.0 5.0 0.333
Result: partially reproduced. Fmr1-null has the larger share of interneurons, as claimed, under every reading: R1, all clusters, 44.5% against 33.5% (z = -2.36, p = 0.018); R2, 48.5% against 35.4% (z = -3.33, p = 0.0009); with the refractory filter p = 0.092 (R1) and 0.023 (R2). The percentages are far from the paper's under every reading: WT 33.0% to 41.4% against 12.4 +/- 1.24%, Fmr1-null 41.0% to 55.7% against 29.3 +/- 2.93%. Interneurons here are clusters above 2 AP/s, without the waveform criterion, among every sorted cluster. Few clusters have refractory violations (70 of 774 have more than 0.5% of intervals below 2 ms), and the filter changes the percentages by at most 4.9 points, so this filter does not remove the excess. It is a weak proxy for the paper's isolation measures: interneuron rates are still 4.34 (WT) and 6.95 (Fmr1-null) AP/s at the median in R1, so whether poorly isolated clusters explain the excess cannot be told from these data.
place-cell-proportion¶
Place cells are less prevalent in Fmr1-null mice (141, 42.2%) than in WT mice (108, 47.8%).
Criterion: counts and percentages within tolerance under one denominator (all cells, or the principal cells), and a two-proportion z test with p < 0.05 (one-sided; two-sided as a variant), Fmr1-null lower.
In [10]:
rows = []
for rname, t in [("R1", R1), ("R2", R2), ("R1 32x32", R1_code), ("R2 32x32", R2_code),
("R1 refractory filter", R1_iso), ("R2 refractory filter", R2_iso)]:
for col in [c for c in ["cls", "cls_rate"] if c in t]:
for den_name, den in [("all clusters", None), ("classified cells", ["interneuron", "place", "nonspatial"]),
("principal cells", ["place", "nonspatial"])]:
if col == "cls_rate" and den == ["interneuron", "place", "nonspatial"]:
continue
o, zz, pp = proportions(t, col, "place", den)
rows.append(dict(reading=rname, classes="rate+burst" if col == "cls" else "rate only",
denominator=den_name, WT=o["WT"][0], WT_n=o["WT"][1], WT_pct=100 * o["WT"][0] / o["WT"][1],
KO=o["Fmr1-null"][0], KO_n=o["Fmr1-null"][1],
KO_pct=100 * o["Fmr1-null"][0] / o["Fmr1-null"][1], z=zz,
p_one_sided=stats.norm.sf(zz), p_two_sided=pp))
pc = pd.DataFrame(rows)
pc["counts_in_tol"] = ((pc.WT - 108).abs() <= 10.8) & ((pc.KO - 141).abs() <= 14.1)
pc["pct_in_tol"] = ((pc.WT_pct - 47.8).abs() <= 4.78) & ((pc.KO_pct - 42.2).abs() <= 4.22)
print(pc.round(4).to_string())
print("\nBy mouse (R1, rate+burst): place cells / all clusters")
print(R1.groupby(["genotype", "mouse"]).apply(lambda d: f"{(d.cls == 'place').sum()}/{len(d)}").to_string())
reading classes denominator WT WT_n WT_pct KO KO_n KO_pct z p_one_sided p_two_sided counts_in_tol pct_in_tol
0 R1 rate+burst all clusters 48 203 23.6453 51 238 21.4286 0.5561 0.2891 0.5782 False False
1 R1 rate+burst classified cells 48 176 27.2727 51 205 24.8780 0.5314 0.2976 0.5952 False False
2 R1 rate+burst principal cells 48 108 44.4444 51 99 51.5152 -1.0173 0.8455 0.3090 False False
3 R1 rate only all clusters 51 203 25.1232 54 238 22.6891 0.5982 0.2749 0.5497 False False
4 R1 rate only principal cells 51 135 37.7778 54 132 40.9091 -0.5237 0.6998 0.6005 False False
5 R2 rate+burst all clusters 55 291 18.9003 70 342 20.4678 -0.4937 0.6892 0.6215 False False
6 R2 rate+burst classified cells 55 249 22.0884 70 298 23.4899 -0.3888 0.6513 0.6974 False False
7 R2 rate+burst principal cells 55 146 37.6712 70 132 53.0303 -2.5707 0.9949 0.0101 False False
8 R2 rate only all clusters 59 291 20.2749 73 342 21.3450 -0.3303 0.6294 0.7412 False False
9 R2 rate only principal cells 59 188 31.3830 73 176 41.4773 -2.0018 0.9773 0.0453 False False
10 R1 32x32 rate+burst all clusters 55 203 27.0936 63 238 26.4706 0.1473 0.4414 0.8829 False False
11 R1 32x32 rate+burst classified cells 55 176 31.2500 63 205 30.7317 0.1091 0.4566 0.9131 False False
12 R1 32x32 rate+burst principal cells 55 108 50.9259 63 99 63.6364 -1.8452 0.9675 0.0650 False False
13 R2 32x32 rate+burst all clusters 65 291 22.3368 81 342 23.6842 -0.4011 0.6558 0.6884 False False
14 R2 32x32 rate+burst classified cells 65 249 26.1044 81 298 27.1812 -0.2835 0.6116 0.7768 False False
15 R2 32x32 rate+burst principal cells 65 146 44.5205 81 132 61.3636 -2.8083 0.9975 0.0050 False False
16 R1 refractory filter rate+burst all clusters 47 194 24.2268 51 222 22.9730 0.3006 0.3818 0.7637 False False
17 R1 refractory filter rate+burst classified cells 47 168 27.9762 51 189 26.9841 0.2097 0.4170 0.8339 False False
18 R1 refractory filter rate+burst principal cells 47 104 45.1923 51 98 52.0408 -0.9734 0.8348 0.3304 False False
19 R1 refractory filter rate only all clusters 50 194 25.7732 54 222 24.3243 0.3405 0.3668 0.7335 False False
20 R1 refractory filter rate only principal cells 50 130 38.4615 54 131 41.2214 -0.4553 0.6756 0.6489 False False
21 R2 refractory filter rate+burst all clusters 54 276 19.5652 70 307 22.8013 -0.9534 0.8298 0.3404 False False
22 R2 refractory filter rate+burst classified cells 54 235 22.9787 70 263 26.6160 -0.9370 0.8256 0.3487 False False
23 R2 refractory filter rate+burst principal cells 54 140 38.5714 70 129 54.2636 -2.5794 0.9951 0.0099 False False
24 R2 refractory filter rate only all clusters 58 276 21.0145 73 307 23.7785 -0.7983 0.7877 0.4247 False False
25 R2 refractory filter rate only principal cells 58 181 32.0442 73 173 42.1965 -1.9776 0.9760 0.0480 False False
By mouse (R1, rate+burst): place cells / all clusters
genotype mouse
Fmr1-null M10 0/5
M18 8/40
M29 8/41
M33 3/20
M34 16/72
M37 14/48
M9 2/12
WT M1 4/9
M11 1/8
M16 14/77
M28 6/24
M3 1/3
M32 13/39
M35 7/28
M36 2/15
Result: not reproduced. No reading gives the paper's counts or percentages: with all clusters as the denominator there are 48 WT (23.6%) and 51 Fmr1-null (21.4%) place cells in R1, and 55 (18.9%) and 70 (20.5%) in R2, against 108 (47.8%) and 141 (42.2%); with principal cells as the denominator WT is 44.4% and Fmr1-null 51.5% (R1). No one-sided test over the 26 readings, including the 32 x 32 grid and the refractory filter, finds Fmr1-null lower (smallest one-sided p 0.27); with principal cells in R2 Fmr1-null is higher (53.0% against 37.7%, two-sided p = 0.010; on the 32 x 32 grid 61.4% against 44.5%, p = 0.005).
figs3-nonspatial-properties-similar¶
Waveform width, burst probability, interspike interval and overall firing rate do not differ between genotypes.
Criterion: a two-way genotype x cell class ANOVA per measure (Figure S3 legend); burst probability is the fraction of ISIs below 10 ms; ISI as the mean per cell (median as a variant); reproduced when no genotype effect has p < 0.05. Waveforms are not in the data set, so width is not tested. Cells: the classified clusters (place, non-spatial pyramidal, interneuron) of R1; R2 shown as a variant.
In [11]:
def anova_rows(t, measures, label):
rows = []
for m in measures:
a = hc16.anova2(t, m)
rows.append(dict(reading=label, measure=m, n=int(t[m].notna().sum()),
F_genotype=a.loc["C(genotype)", "F"], df_resid=int(a.loc["Residual", "df"]),
p_genotype=a.loc["C(genotype)", "PR(>F)"], eta_p2_genotype=a.loc["C(genotype)", "eta_p2"],
p_class=a.loc["C(cls)", "PR(>F)"], p_interaction=a.loc["C(genotype):C(cls)", "PR(>F)"]))
return rows
rows = []
for label, t in [("R1", R1), ("R2", R2), ("R1 refractory filter", R1_iso), ("R2 refractory filter", R2_iso)]:
tc = t[t.cls.isin(["place", "nonspatial", "interneuron"])].copy()
rows += anova_rows(tc, ["burst_prob", "isi_mean", "isi_median", "rate"], label)
f3 = pd.DataFrame(rows)
print(f3.round(4).to_string())
tc = R1[R1.cls.isin(["place", "nonspatial", "interneuron"])]
print("\nR1 means by genotype and class:")
print(tc.groupby(["cls", "genotype"])[["burst_prob", "isi_mean", "isi_median", "rate"]].agg(["mean", "count"]).round(3).to_string())
print("\nR1 by mouse: mean burst probability and rate of classified clusters")
print(tc.groupby(["genotype", "mouse"])[["burst_prob", "rate"]].mean().round(3).to_string())
reading measure n F_genotype df_resid p_genotype eta_p2_genotype p_class p_interaction
0 R1 burst_prob 381 2.3626 375 0.1251 0.0063 0.0 0.0775
1 R1 isi_mean 381 0.7112 375 0.3996 0.0019 0.0 0.4041
2 R1 isi_median 381 0.2353 375 0.6279 0.0006 0.0 0.7448
3 R1 rate 381 8.2848 375 0.0042 0.0216 0.0 0.0086
4 R2 burst_prob 547 4.8146 541 0.0286 0.0088 0.0 0.0989
5 R2 isi_mean 547 0.7055 541 0.4013 0.0013 0.0 0.4043
6 R2 isi_median 547 0.3706 541 0.5429 0.0007 0.0 0.5839
7 R2 rate 547 12.1301 541 0.0005 0.0219 0.0 0.0024
8 R1 refractory filter burst_prob 357 3.1192 351 0.0782 0.0088 0.0 0.0112
9 R1 refractory filter isi_mean 357 0.7669 351 0.3818 0.0022 0.0 0.3834
10 R1 refractory filter isi_median 357 0.2954 351 0.5871 0.0008 0.0 0.6958
11 R1 refractory filter rate 357 3.3099 351 0.0697 0.0093 0.0 0.1460
12 R2 refractory filter burst_prob 498 8.5237 492 0.0037 0.0170 0.0 0.0061
13 R2 refractory filter isi_mean 498 0.4894 492 0.4845 0.0010 0.0 0.5367
14 R2 refractory filter isi_median 498 0.3585 492 0.5496 0.0007 0.0 0.5991
15 R2 refractory filter rate 498 3.6393 492 0.0570 0.0073 0.0 0.1564
R1 means by genotype and class:
burst_prob isi_mean isi_median rate
mean count mean count mean count mean count
cls genotype
interneuron Fmr1-null 0.148 106 0.179 106 0.075 106 11.473 106
WT 0.185 68 0.239 68 0.072 68 7.270 68
nonspatial Fmr1-null 0.153 48 14.297 48 4.908 48 0.531 48
WT 0.131 60 20.313 60 6.373 60 0.373 60
place Fmr1-null 0.200 51 1.847 51 0.161 51 0.915 51
WT 0.224 48 1.681 48 0.160 48 0.956 48
R1 by mouse: mean burst probability and rate of classified clusters
burst_prob rate
genotype mouse
Fmr1-null M10 0.242 1.735
M18 0.198 3.812
M29 0.099 7.586
M33 0.149 2.320
M34 0.152 7.140
M37 0.188 7.741
M9 0.185 5.104
WT M1 0.146 1.940
M11 0.222 2.759
M16 0.166 3.553
M28 0.132 3.240
M3 0.112 2.278
M32 0.245 1.941
M35 0.170 4.203
M36 0.136 4.224
Result: partially reproduced. In R1 (381 classified clusters) burst probability and the ISI show no genotype effect (burst p = 0.13; mean ISI p = 0.40; median ISI p = 0.63), but the rate does (F = 8.28, p = 0.0042, partial eta squared 0.022, with an interaction p = 0.0086; Fmr1-null interneurons fire at 11.5 against 7.3 AP/s). In R2 burst probability also differs (p = 0.029) and the rate (p = 0.0005). With the refractory filter the rate effect falls below significance (R1 p = 0.070, R2 p = 0.057), so that in R1 no measure differs, but in R2 burst probability does (p = 0.0037). Waveform width is not in the data set and is not tested.
fig2-spatial-quality-similar¶
In a fixed environment, spatial firing quality (coherence, proportion of active pixels, information content, number of place fields, in-field rate) differs between place cells and interneurons but not between genotypes.
Criterion: a two-way genotype x cell class ANOVA per measure on cells recorded in the box; the box is open (large, small or both), and only the large box is in the data set. Fields are 4-connected pixels above (a) the cell's mean rate or (b) 20% of its peak rate, at least 9 pixels. Reproduced when no genotype effect has p < 0.05 and the class effect has p < 0.05 for coherence, active pixels, information and in-field rate.
In [12]:
box = cells[(cells.env == "big_box") & cells.cls.isin(["place", "nonspatial", "interneuron"])].copy()
print("Large box cells by genotype and class:")
print(box.groupby(["genotype", "cls"]).size().unstack())
print("Mice:", box.groupby("genotype").mouse.nunique().to_dict())
measures = ["coherence", "active", "info", "nfields_mean", "infield_mean", "nfields_20pk", "infield_20pk"]
f2 = pd.DataFrame(anova_rows(box, measures, "large box"))
box_code = cells_code[(cells_code.env == "big_box") & cells_code.cls.isin(["place", "nonspatial", "interneuron"])].copy()
box_iso = box[box.refractory <= ISO]
f2 = pd.concat([f2, pd.DataFrame(anova_rows(box_code, measures[:-1], "32x32 grid")),
pd.DataFrame(anova_rows(box_iso, measures, "refractory filter"))], ignore_index=True)
print(f2.round(4).to_string())
print("32x32 grid, large box cells by genotype and class:")
print(box_code.groupby(["genotype", "cls"]).size().unstack().to_string())
print("\nMeans by class and genotype:")
print(box.groupby(["cls", "genotype"])[measures].mean().round(3).to_string())
print("\nRanges:", box[measures].describe().loc[["min", "max"]].round(3).to_string())
fig, axs = plt.subplots(1, 5, figsize=(15, 3))
for ax, m in zip(axs, ["coherence", "active", "info", "nfields_mean", "infield_mean"]):
for i, (c, g) in enumerate([(c, g) for c in ["place", "nonspatial", "interneuron"] for g in ["WT", "Fmr1-null"]]):
v = box[(box.cls == c) & (box.genotype == g)][m].dropna()
ax.errorbar(i, v.mean(), v.sem(), fmt="o", color="k" if g == "WT" else "r")
ax.set_xticks(range(6)); ax.set_xticklabels(["PC WT", "PC KO", "NS WT", "NS KO", "IN WT", "IN KO"], rotation=60)
ax.set_title(m)
plt.tight_layout(); plt.show()
Large box cells by genotype and class:
cls interneuron nonspatial place
genotype
Fmr1-null 46 35 36
WT 33 29 34
Mice: {'Fmr1-null': 7, 'WT': 8}
reading measure n F_genotype df_resid p_genotype eta_p2_genotype p_class p_interaction
0 large box coherence 213 9.0099 207 0.0030 0.0417 0.0000 0.0077
1 large box active 213 18.7084 207 0.0000 0.0829 0.0000 0.0098
2 large box info 213 7.2717 207 0.0076 0.0339 0.0000 0.1789
3 large box nfields_mean 213 6.6162 207 0.0108 0.0310 0.0000 0.0951
4 large box infield_mean 199 2.5431 193 0.1124 0.0130 0.0000 0.0541
5 large box nfields_20pk 213 0.4994 207 0.4805 0.0024 0.0000 0.6933
6 large box infield_20pk 100 0.6435 94 0.4245 0.0068 0.0004 0.5706
7 32x32 grid coherence 213 13.7992 207 0.0003 0.0625 0.0000 0.0036
8 32x32 grid active 213 10.2814 207 0.0016 0.0473 0.0000 0.0998
9 32x32 grid info 213 7.3377 207 0.0073 0.0342 0.0000 0.0294
10 32x32 grid nfields_mean 213 3.4395 207 0.0651 0.0163 0.0000 0.1139
11 32x32 grid infield_mean 200 2.5852 194 0.1095 0.0132 0.0000 0.0558
12 32x32 grid nfields_20pk 213 0.0778 207 0.7805 0.0004 0.0000 0.8336
13 refractory filter coherence 205 7.0948 199 0.0084 0.0344 0.0000 0.0292
14 refractory filter active 205 15.4738 199 0.0001 0.0721 0.0000 0.0258
15 refractory filter info 205 6.2358 199 0.0133 0.0304 0.0000 0.1783
16 refractory filter nfields_mean 205 4.7762 199 0.0300 0.0234 0.0000 0.2002
17 refractory filter infield_mean 192 1.6125 186 0.2057 0.0086 0.0000 0.1314
18 refractory filter nfields_20pk 205 0.0444 199 0.8333 0.0002 0.0000 0.9642
19 refractory filter infield_20pk 96 0.4618 90 0.4985 0.0051 0.0008 0.6211
32x32 grid, large box cells by genotype and class:
cls interneuron nonspatial place
genotype
Fmr1-null 46 25 46
WT 33 23 40
Means by class and genotype:
coherence active info nfields_mean infield_mean nfields_20pk infield_20pk
cls genotype
interneuron Fmr1-null 0.416 0.821 0.402 8.913 14.963 3.261 19.288
WT 0.615 0.647 0.819 6.242 10.609 2.727 16.444
nonspatial Fmr1-null 0.201 0.222 2.501 3.686 2.135 0.029 5.920
WT 0.206 0.152 3.214 3.000 2.113 0.034 7.140
place Fmr1-null 0.623 0.261 2.266 3.083 3.948 0.861 11.307
WT 0.651 0.240 2.315 2.824 4.531 0.824 11.900
Ranges: coherence active info nfields_mean infield_mean nfields_20pk infield_20pk
min -0.055 0.003 0.041 0.0 0.493 0.0 3.720
max 1.265 0.996 8.740 20.0 54.393 14.0 53.765
Result: not reproduced. In the large box (213 classified clusters) the genotype effect is significant for coherence (F = 9.01, p = 0.0030), active pixels (F = 18.7, p < 0.0001), information (F = 7.27, p = 0.0076) and the number of fields above the mean rate (p = 0.011), with interactions for coherence (p = 0.0077) and active pixels (p = 0.0098); the differences lie mostly in interneurons (coherence 0.42 Fmr1-null against 0.62 WT). In-field rate (p = 0.11 and 0.42 under the two field readings) and the count of fields above 20% of the peak (p = 0.48) show none. The class effect is significant for every measure (p <= 0.0004). The genotype effects on coherence, active pixels and information remain on the 32 x 32 grid (p = 0.0003, 0.0016, 0.0073) and with the refractory filter (p = 0.0084, 0.0001, 0.013).
tables1-place-cell-properties-similar¶
Firing rate, in-field rate, spatial coherence, information content and proportion of active pixels of place cells, and the firing rate of interneurons, do not differ between WT and Fmr1-null mice.
Criterion: a t test per measure, or Mann-Whitney when a sample fails Shapiro-Wilk (Table S1 legend); reproduced when no comparison has p < 0.05. Cells and box as in fig2 (large box), in-field rate under both field readings. Means are compared with Table S1 (tolerances its SEMs).
In [13]:
expected = {("place", "rate"): (1.49, 0.15, 1.27, 0.12), ("place", "infield_mean"): (2.61, 0.2, 2.2, 0.13),
("place", "infield_20pk"): (2.61, 0.2, 2.2, 0.13), ("place", "coherence"): (0.61, 0.02, 0.57, 0.02),
("place", "info"): (1.57, 0.09, 1.5, 0.09), ("place", "active"): (0.45, 0.03, 0.46, 0.03),
("interneuron", "rate"): (16.61, 2.8, 16.02, 1.42)}
rows = []
for variant, bx in [("main", box), ("32x32 grid", box_code), ("refractory filter", box_iso)]:
for (c, m), (ew, sw, ek, sk) in expected.items():
a = bx[(bx.cls == c) & (bx.genotype == "WT")][m]
b = bx[(bx.cls == c) & (bx.genotype == "Fmr1-null")][m]
if min(a.notna().sum(), b.notna().sum()) < 3:
continue
r = hc16.two_group(a, b)
rows.append(dict(variant=variant, cls=c, measure=m, n_WT=r["n_a"], WT=r["mean_a"], WT_sem=r["sem_a"], n_KO=r["n_b"], KO=r["mean_b"],
KO_sem=r["sem_b"], t=r["t"], p_t=r["p_t"], U=r["U"], p_U=r["p_U"], test=r["test_table"],
p_table=r["p_table"], d=r["d"], WT_in_tol=abs(r["mean_a"] - ew) <= sw, KO_in_tol=abs(r["mean_b"] - ek) <= sk))
s1 = pd.DataFrame(rows)
print(s1.round(4).to_string())
print("\nPlace cells by mouse (large box): n, mean rate, coherence, information")
print(box[box.cls == "place"].groupby(["genotype", "mouse"])[["rate", "coherence", "info"]].agg(["count", "mean"]).round(3).to_string())
variant cls measure n_WT WT WT_sem n_KO KO KO_sem t p_t U p_U test p_table d WT_in_tol KO_in_tol
0 main place rate 34 0.9740 0.0924 36 0.8874 0.0824 0.7016 0.4853 659.0 0.5848 t 0.4853 0.1678 False False
1 main place infield_mean 34 4.5313 0.3090 36 3.9482 0.2549 1.4623 0.1483 722.0 0.1982 t 0.1483 0.3497 False False
2 main place infield_20pk 21 11.8996 1.1021 21 11.3072 0.9762 0.4024 0.6896 230.0 0.8209 t 0.6896 0.1242 False False
3 main place coherence 34 0.6515 0.0297 36 0.6226 0.0289 0.6963 0.4886 685.0 0.3942 U 0.3942 0.1665 False False
4 main place info 34 2.3154 0.1118 36 2.2662 0.1520 0.2586 0.7967 674.0 0.4699 U 0.4699 0.0618 False False
5 main place active 34 0.2398 0.0187 36 0.2613 0.0202 -0.7797 0.4383 538.0 0.3878 t 0.4383 -0.1865 False False
6 main interneuron rate 33 5.3984 0.9901 46 10.3430 1.5281 -2.4827 0.0152 501.0 0.0105 U 0.0105 -0.5664 False False
7 32x32 grid place rate 40 0.9345 0.0840 46 0.8324 0.0740 0.9160 0.3623 1015.0 0.4132 U 0.4132 0.1980 False False
8 32x32 grid place infield_mean 40 3.3896 0.2666 46 2.8773 0.2029 1.5507 0.1247 1069.0 0.1985 t 0.1247 0.3352 False False
9 32x32 grid place infield_20pk 33 6.3145 0.4375 33 5.4814 0.5220 1.2233 0.2257 656.0 0.1546 U 0.1546 0.3012 False False
10 32x32 grid place coherence 40 0.8650 0.0415 46 0.7940 0.0384 1.2550 0.2130 1074.0 0.1838 U 0.1838 0.2713 False False
11 32x32 grid place info 40 1.7444 0.0913 46 1.7356 0.1254 0.0554 0.9560 988.0 0.5589 U 0.5589 0.0120 False False
12 32x32 grid place active 40 0.3990 0.0255 46 0.4136 0.0255 -0.4040 0.6872 866.0 0.6432 U 0.6432 -0.0873 False False
13 32x32 grid interneuron rate 33 5.3984 0.9901 46 10.3430 1.5281 -2.4827 0.0152 501.0 0.0105 U 0.0105 -0.5664 False False
14 refractory filter place rate 33 0.9850 0.0945 36 0.8874 0.0824 0.7811 0.4375 646.0 0.5362 t 0.4375 0.1882 False False
15 refractory filter place infield_mean 33 4.5377 0.3184 36 3.9482 0.2549 1.4563 0.1500 700.0 0.2050 t 0.1500 0.3510 False False
16 refractory filter place infield_20pk 21 11.8996 1.1021 21 11.3072 0.9762 0.4024 0.6896 230.0 0.8209 t 0.6896 0.1242 False False
17 refractory filter place coherence 33 0.6559 0.0303 36 0.6226 0.0289 0.7952 0.4293 673.0 0.3457 U 0.3457 0.1916 False False
18 refractory filter place info 33 2.3010 0.1142 36 2.2662 0.1520 0.1804 0.8574 647.0 0.5283 U 0.5283 0.0435 False False
19 refractory filter place active 33 0.2436 0.0188 36 0.2613 0.0202 -0.6376 0.5259 533.0 0.4674 t 0.5259 -0.1537 False False
20 refractory filter interneuron rate 32 5.4539 1.0200 42 9.6665 1.6050 -2.0601 0.0430 479.0 0.0357 U 0.0357 -0.4834 False False
Place cells by mouse (large box): n, mean rate, coherence, information
rate coherence info
count mean count mean count mean
genotype mouse
Fmr1-null M18 4 0.410 4 0.581 4 3.259
M29 7 1.127 7 0.581 7 1.804
M33 3 0.689 3 0.624 3 2.249
M34 9 1.047 9 0.729 9 2.210
M37 11 0.816 11 0.587 11 2.363
M9 2 0.974 2 0.564 2 1.643
WT M1 4 1.065 4 0.596 4 2.045
M11 1 0.455 1 0.636 1 2.647
M16 10 0.938 10 0.718 10 2.670
M28 4 0.988 4 0.632 4 2.400
M3 1 0.459 1 0.810 1 1.825
M32 9 1.155 9 0.566 9 2.064
M35 5 0.843 5 0.704 5 2.237
Result: partially reproduced. Among large-box place cells (34 WT, 36 Fmr1-null) no measure differs: firing rate p = 0.49, in-field rate p = 0.15 (mean-rate fields) and 0.69 (20% of peak), coherence p = 0.39 (Mann-Whitney), information p = 0.47 (Mann-Whitney), active pixels p = 0.44. The interneuron rates differ (WT 5.40, Fmr1-null 10.34 AP/s; Mann-Whitney p = 0.0105; with the refractory filter p = 0.036). No mean is within Table S1's SEM: information is 2.32 and 2.27 bits/spike against 1.57 and 1.5, active pixels 0.24 and 0.26 against 0.45 and 0.46, firing rate 0.97 and 0.89 against 1.49 and 1.27 AP/s. On the 32 x 32 grid (40 and 46 place cells) no place cell measure differs either (p = 0.12 to 0.64), and information (1.74 and 1.74 bits/spike) and active pixels (0.40 and 0.41) come nearer Table S1 without reaching it.
tables1-nonspatial-active-pixels¶
Non-spatial pyramidal cells of Fmr1-null mice discharge in a larger proportion of the environment than those of WT mice (0.40 against 0.26).
Criterion: both values within tolerance (0.26 +/- 0.03, 0.40 +/- 0.05) and the t test (or Mann-Whitney, by Shapiro-Wilk) p < 0.05 with Fmr1-null higher. Large box cells.
In [14]:
a = box[(box.cls == "nonspatial") & (box.genotype == "WT")]["active"]
b = box[(box.cls == "nonspatial") & (box.genotype == "Fmr1-null")]["active"]
r = hc16.two_group(a, b)
print({k: (round(v, 4) if isinstance(v, float) else v) for k, v in r.items()})
print("WT within 0.26 +/- 0.03:", abs(r["mean_a"] - 0.26) <= 0.03, "; Fmr1-null within 0.40 +/- 0.05:", abs(r["mean_b"] - 0.40) <= 0.05)
for label, t in [("R1 (all first stationary sessions)", R1), ("R2", R2), ("large box, 32x32 grid", box_code),
("large box, refractory filter", box_iso)]:
rr = hc16.two_group(t[(t.cls == "nonspatial") & (t.genotype == "WT")]["active"],
t[(t.cls == "nonspatial") & (t.genotype == "Fmr1-null")]["active"])
print(f"Variant {label}: WT {rr['mean_a']:.3f} (n={rr['n_a']}), Fmr1-null {rr['mean_b']:.3f} (n={rr['n_b']}), "
f"t p = {rr['p_t']:.4f}, U p = {rr['p_U']:.4f}, d = {rr['d']:.2f}")
print("By mouse (large box):")
print(box[box.cls == "nonspatial"].groupby(["genotype", "mouse"]).active.agg(["count", "mean"]).round(3).to_string())
{'n_a': 29, 'n_b': 35, 'mean_a': np.float64(0.1519), 'sem_a': np.float64(0.0275), 'mean_b': np.float64(0.2218), 'sem_b': np.float64(0.0309), 't': np.float64(-1.6582), 'df': 62, 'p_t': np.float64(0.1023), 'U': np.float64(385.0), 'p_U': np.float64(0.0999), 'shapiro_a': np.float64(0.0014), 'shapiro_b': np.float64(0.007), 'p_table': np.float64(0.0999), 'test_table': 'U', 'd': -0.4164}
WT within 0.26 +/- 0.03: False ; Fmr1-null within 0.40 +/- 0.05: False
Variant R1 (all first stationary sessions): WT 0.200 (n=60), Fmr1-null 0.250 (n=48), t p = 0.2255, U p = 0.0944, d = -0.24
Variant R2: WT 0.228 (n=91), Fmr1-null 0.280 (n=62), t p = 0.1529, U p = 0.0986, d = -0.24
Variant large box, 32x32 grid: WT 0.241 (n=23), Fmr1-null 0.355 (n=25), t p = 0.1135, U p = 0.0987, d = -0.47
Variant large box, refractory filter: WT 0.147 (n=27), Fmr1-null 0.222 (n=35), t p = 0.0861, U p = 0.0911, d = -0.45
By mouse (large box):
count mean
genotype mouse
Fmr1-null M10 2 0.316
M18 5 0.206
M29 2 0.287
M33 7 0.145
M34 8 0.273
M37 8 0.166
M9 3 0.334
WT M1 1 0.205
M11 2 0.071
M16 12 0.115
M28 2 0.470
M3 1 0.009
M32 6 0.161
M35 5 0.152
Result: not reproduced. Large-box non-spatial pyramidal clusters discharge in 0.152 (WT, n = 29) and 0.222 (Fmr1-null, n = 35) of the visited pixels. Fmr1-null is higher, but both samples fail Shapiro-Wilk and Mann-Whitney gives p = 0.10 (t test p = 0.10), and neither value is within tolerance of 0.26 and 0.40. R1 and R2 give p = 0.09 to 0.23. On the 32 x 32 grid both values are within tolerance (0.241 and 0.355, 23 and 25 cells) but the difference is not significant (Mann-Whitney p = 0.099); with the refractory filter p = 0.091.
Passes and overdispersion¶
STAR Methods, "Overdispersion": for each 5 s interval, exp is the sum over position samples of the time spent times the session's time-averaged rate at that pixel, z = (obs - exp)/sqrt(exp), and the interval is a pass when exp is at least the cell's average rate. The claims file reads this literally, a count against a rate in AP/s; the count expected at the average rate over the 5 s (exp at least the average rate x 5 s) is a variant, and which is meant is asked of a reader (questions/crcns-hc-16.yaml). Intervals are consecutive from the session start (the claims file's reading); a 1 s sliding step is a variant. Place cells of the stationary sessions (large box, stationary arena 1 and 2), each session on its own (R2); R1 (each cell once, from its first session) is a variant.
In [15]:
def pass_table(t, step_s=None, rule="rate", index_fn=hc16.pixel_index):
out = {}
for sid, g in t[t.cls == "place"].groupby("session"):
s = S[sid]
ish = index_fn(s, "Rm")
for c in g.cell:
out[(sid, c)] = hc16.passes(s, c, step_s=step_s, idx_shape=ish, rule=rule)
return out
P_R2 = pass_table(R2)
P_R2_step1 = pass_table(R2, step_s=1.0)
P_R2_count = pass_table(R2, rule="count")
pz = pd.DataFrame([dict(session=k[0], cell=k[1], genotype=S[k[0]].genotype, mouse=S[k[0]].mouse, z=zz)
for k, v in P_R2.items() for zz in v["z"][v["is_pass"]]])
pz["in_R1"] = pz.session.isin(set(R1.session))
pz_count = pd.DataFrame([dict(session=k[0], cell=k[1], genotype=S[k[0]].genotype, mouse=S[k[0]].mouse, z=zz)
for k, v in P_R2_count.items() for zz in v["z"][v["is_pass"]]])
pz_count["in_R1"] = pz_count.session.isin(set(R1.session))
for name, P in [("literal (exp >= mean rate)", P_R2), ("count (exp >= mean rate x 5 s)", P_R2_count)]:
f = np.array([v["is_pass"].mean() for v in P.values()]); n = np.array([v["is_pass"].sum() for v in P.values()])
print(f"Pass rule {name}: fraction of 5 s intervals that pass, median {np.median(f):.3f} (range {f.min():.3f} to {f.max():.3f}); "
f"passes per cell median {np.median(n):.0f} (range {n.min()} to {n.max()})")
print("Place cells with passes:", pz.groupby("genotype")[["session", "cell"]].apply(lambda d: len(d.drop_duplicates())).to_dict())
print("z over passes: range %.2f to %.2f" % (pz.z.min(), pz.z.max()))
Pass rule literal (exp >= mean rate): fraction of 5 s intervals that pass, median 0.763 (range 0.329 to 1.000); passes per cell median 241 (range 93 to 359)
Pass rule count (exp >= mean rate x 5 s): fraction of 5 s intervals that pass, median 0.331 (range 0.167 to 0.460); passes per cell median 109 (range 49 to 165)
Place cells with passes: {'Fmr1-null': 70, 'WT': 55}
z over passes: range -8.48 to 16.87
overdispersion-passes¶
Overdispersion is lower (spatial firing more reliable) in Fmr1-null mice when the unit is a pass through a place field (variance 4.89 against 5.81).
Criterion: both variances within tolerance (5.81 +/- 0.58, 4.89 +/- 0.49), and the F test and Levene's test p < 0.05 with Fmr1-null lower. Pass counts expected 6,776 and 9,102 (+/- 10%).
In [16]:
def od_passes(zz, label):
a = zz[zz.genotype == "WT"].z.to_numpy(); b = zz[zz.genotype == "Fmr1-null"].z.to_numpy()
F = a.var(ddof=1) / b.var(ddof=1)
p_one = stats.f.sf(F, len(a) - 1, len(b) - 1)
lev = stats.levene(a, b, center="mean")
return dict(reading=label, n_WT=len(a), var_WT=a.var(ddof=1), n_KO=len(b), var_KO=b.var(ddof=1), F=F,
p_F_one_sided=p_one, p_F_two_sided=2 * min(p_one, 1 - p_one), levene_F=lev.statistic, levene_p=lev.pvalue,
WT_in_tol=abs(a.var(ddof=1) - 5.81) <= 0.58, KO_in_tol=abs(b.var(ddof=1) - 4.89) <= 0.49,
passes_in_tol=(abs(len(a) - 6776) <= 678) and (abs(len(b) - 9102) <= 910))
pz1 = pd.DataFrame([dict(session=k[0], genotype=S[k[0]].genotype, z=zz)
for k, v in P_R2_step1.items() for zz in v["z"][v["is_pass"]]])
odp = pd.DataFrame([od_passes(pz, "R2, consecutive 5 s"), od_passes(pz[pz.in_R1], "R1, consecutive 5 s"),
od_passes(pz1, "R2, 5 s every 1 s"),
od_passes(pz_count, "R2, count rule"), od_passes(pz_count[pz_count.in_R1], "R1, count rule")])
print(odp.round(4).to_string())
print("\nBy mouse (R2): passes and variance of z")
print(pz.groupby(["genotype", "mouse"]).z.agg(["count", "var"]).round(3).to_string())
print("\nBy mouse (R2, count rule): passes and variance of z")
print(pz_count.groupby(["genotype", "mouse"]).z.agg(["count", "var"]).round(3).to_string())
reading n_WT var_WT n_KO var_KO F p_F_one_sided p_F_two_sided levene_F levene_p WT_in_tol KO_in_tol passes_in_tol
0 R2, consecutive 5 s 12802 3.6903 18521 3.5493 1.0397 0.0082 0.0164 0.0000 0.9977 False False False
1 R1, consecutive 5 s 10747 3.5339 12815 3.8267 0.9235 1.0000 0.0000 20.8587 0.0000 False False False
2 R2, 5 s every 1 s 63865 3.6561 92254 3.5918 1.0179 0.0074 0.0148 0.2564 0.6126 False False False
3 R2, count rule 5442 5.2789 8250 4.6635 1.1319 0.0000 0.0000 11.6845 0.0006 True True False
4 R1, count rule 4609 4.9971 5643 5.1627 0.9679 0.8769 0.2462 4.0172 0.0451 False True False
By mouse (R2): passes and variance of z
count var
genotype mouse
Fmr1-null M18 3452 4.553
M29 2095 3.697
M33 838 2.340
M34 7116 2.971
M37 4494 3.968
M9 526 2.560
WT M1 971 4.099
M11 258 3.604
M16 3345 3.179
M28 1422 5.027
M3 278 3.140
M32 3562 3.808
M35 2262 3.177
M36 704 4.168
By mouse (R2, count rule): passes and variance of z
count var
genotype mouse
Fmr1-null M18 1502 5.825
M29 895 5.639
M33 334 3.412
M34 3293 3.738
M37 1940 5.384
M9 286 2.763
WT M1 425 5.396
M11 124 4.180
M16 1503 4.180
M28 546 8.717
M3 121 4.580
M32 1511 5.591
M35 954 4.192
M36 258 7.293
Result: not reproduced, on the claims file's pass rule. In R2 there are 12,802 WT and 18,521 Fmr1-null passes (expected 6,776 and 9,102) with variances of z 3.69 and 3.55 (expected 5.81 and 4.89). The F test finds the WT variance larger (F = 1.04, one-sided p = 0.008) but Levene's test does not (p = 0.998). In R1 the Fmr1-null variance is the larger (3.83 against 3.53). With 5 s intervals every 1 s the result is as in R2 (Levene p = 0.61). The literal rule passes a median 76% of a cell's 5 s intervals. With the count rule (33% of intervals pass) R2 gives 5,442 and 8,250 passes, variances 5.28 and 4.66, both within tolerance, and Fmr1-null lower by the F test (F = 1.13, one-sided p < 0.0001) and Levene's test (p = 0.0006), which meets the criterion; the WT pass count is outside 6,776 +/- 678. In R1 the count rule gives Fmr1-null the larger variance (5.16 against 5.00). Which rule the paper means is asked of a reader.
overdispersion-cells¶
Overdispersion is marginally lower in Fmr1-null mice when the unit is a place cell (5.02 against 5.98).
Criterion: per cell, the variance of z over its passes; mean +/- SEM per genotype; minimum passes per cell 2, 10 or 30; reproduced when under one reading both means are within tolerance (5.98 +/- 0.44, 5.02 +/- 0.3) and the one-tailed t test gives p < 0.05, Fmr1-null lower.
In [17]:
pc_var = pz.groupby(["session", "cell", "genotype", "in_R1"]).z.agg(["count", "var"]).reset_index()
pc_var_count = pz_count.groupby(["session", "cell", "genotype", "in_R1"]).z.agg(["count", "var"]).reset_index()
rows = []
for label, d0 in [("R2", pc_var), ("R1", pc_var[pc_var.in_R1]), ("R2, count rule", pc_var_count),
("R1, count rule", pc_var_count[pc_var_count.in_R1])]:
for mn in [2, 10, 30]:
d = d0[d0["count"] >= mn]
r = hc16.two_group(d[d.genotype == "WT"]["var"], d[d.genotype == "Fmr1-null"]["var"], alternative="greater")
rows.append(dict(reading=label, min_passes=mn, n_WT=r["n_a"], WT=r["mean_a"], WT_sem=r["sem_a"], n_KO=r["n_b"],
KO=r["mean_b"], KO_sem=r["sem_b"], t=r["t"], p_t_one=r["p_t"], U=r["U"], p_U_one=r["p_U"], d=r["d"],
WT_in_tol=abs(r["mean_a"] - 5.98) <= 0.44, KO_in_tol=abs(r["mean_b"] - 5.02) <= 0.3))
odc = pd.DataFrame(rows)
print(odc.round(4).to_string())
print("Per-cell variance range: %.2f to %.2f" % (pc_var["var"].min(), pc_var["var"].max()))
print("Passes per cell, count rule: min %d, median %.0f" % (pc_var_count["count"].min(), pc_var_count["count"].median()))
reading min_passes n_WT WT WT_sem n_KO KO KO_sem t p_t_one U p_U_one d WT_in_tol KO_in_tol 0 R2 2 55 3.6828 0.2409 70 3.5280 0.2134 0.4810 0.3157 2047.0 0.2728 0.0867 False False 1 R2 10 55 3.6828 0.2409 70 3.5280 0.2134 0.4810 0.3157 2047.0 0.2728 0.0867 False False 2 R2 30 55 3.6828 0.2409 70 3.5280 0.2134 0.4810 0.3157 2047.0 0.2728 0.0867 False False 3 R1 2 48 3.5634 0.2524 51 3.7174 0.2533 -0.4303 0.6660 1163.0 0.6666 -0.0865 False False 4 R1 10 48 3.5634 0.2524 51 3.7174 0.2533 -0.4303 0.6660 1163.0 0.6666 -0.0865 False False 5 R1 30 48 3.5634 0.2524 51 3.7174 0.2533 -0.4303 0.6660 1163.0 0.6666 -0.0865 False False 6 R2, count rule 2 55 5.3285 0.4593 70 4.6888 0.3268 1.1650 0.1231 2125.0 0.1605 0.2099 False False 7 R2, count rule 10 55 5.3285 0.4593 70 4.6888 0.3268 1.1650 0.1231 2125.0 0.1605 0.2099 False False 8 R2, count rule 30 55 5.3285 0.4593 70 4.6888 0.3268 1.1650 0.1231 2125.0 0.1605 0.2099 False False 9 R1, count rule 2 48 5.0217 0.4424 51 5.0779 0.3982 -0.0946 0.5376 1188.0 0.6009 -0.0190 False True 10 R1, count rule 10 48 5.0217 0.4424 51 5.0779 0.3982 -0.0946 0.5376 1188.0 0.6009 -0.0190 False True 11 R1, count rule 30 48 5.0217 0.4424 51 5.0779 0.3982 -0.0946 0.5376 1188.0 0.6009 -0.0190 False True Per-cell variance range: 1.05 to 10.93 Passes per cell, count rule: min 49, median 109
Result: not reproduced. Every place cell has at least 30 passes, so the three minimums give the same sample. R2: WT 3.68 +/- 0.24 (n = 55), Fmr1-null 3.53 +/- 0.21 (n = 70), one-tailed t p = 0.32, Mann-Whitney p = 0.27; R1: 3.56 against 3.72 (p = 0.67). Neither mean is within tolerance of 5.98 and 5.02. With the count rule (at least 49 passes per cell, again one sample) R2 gives 5.33 +/- 0.46 against 4.69 +/- 0.33 (one-tailed t p = 0.12) and R1 5.02 against 5.08 (p = 0.54); in no reading are both means within tolerance.
fig4b-covariance-30-passes and fig4b-covariance-40-passes¶
Pairs of Fmr1-null place cells with overlapping fields covary less in their pass-by-pass discharge than WT pairs (average r 0.04 against 0.09), a difference that disappears after shuffling. / With at least 40 passes, the difference holds and is larger (0.06 against 0.19).
Criterion: for each pair of simultaneously recorded place cells, the Pearson correlation of the two cells' z over the intervals that are passes for both, pairs with at least 30 (40) such passes; t test on Fisher z; averages within tolerance (WT 0.09 +/- 0.03, Fmr1-null 0.04 +/- 0.01; for 40, 0.19 +/- 0.04 and 0.06 +/- 0.01); the shuffle permutes z within each cell of a pair (1 permutation per pair; 100 as a variant).
In [18]:
rng_cov = np.random.default_rng(404)
def pair_cov(P, nshuf=100, rng=rng_cov):
rows = []
by_sess = {}
for (sid, c), v in P.items():
by_sess.setdefault(sid, []).append((c, v))
for sid, lst in by_sess.items():
for i in range(len(lst)):
for j in range(i + 1, len(lst)):
(ca, va), (cb, vb) = lst[i], lst[j]
both = va["is_pass"] & vb["is_pass"]
n = int(both.sum())
if n < 3:
continue
za, zb = va["z"][both], vb["z"][both]
r = stats.pearsonr(za, zb).statistic if za.std() > 0 and zb.std() > 0 else np.nan
sh = [stats.pearsonr(rng.permutation(za), rng.permutation(zb)).statistic for _ in range(nshuf)]
rows.append(dict(session=sid, genotype=S[sid].genotype, mouse=S[sid].mouse, a=ca, b=cb, n=n, r=r,
r_shuf1=sh[0], r_shuf100=np.mean(hc16.fisher(np.array(sh))), r_shuf_all=np.array(sh)))
return pd.DataFrame(rows)
cov = pair_cov(P_R2)
cov["in_R1"] = cov.session.isin(set(R1.session))
cov_count = pair_cov(P_R2_count)
cov_count["in_R1"] = cov_count.session.isin(set(R1.session))
for name, d in [("literal rule", cov), ("count rule", cov_count)]:
print(f"Joint passes per pair, {name}: percentiles 0/10/50/90/100 =", np.percentile(d.n, [0, 10, 50, 90, 100]).round(0).tolist(),
"; pairs with >= 30:", int((d.n >= 30).sum()), "and >= 40:", int((d.n >= 40).sum()), "of", len(d))
print("Pairs with >= 3 joint passes:", cov.groupby("genotype").size().to_dict(), "; joint passes range", cov.n.min(), "to", cov.n.max())
rows = []
for label, d0 in [("R2", cov), ("R1", cov[cov.in_R1]), ("R2, count rule", cov_count), ("R1, count rule", cov_count[cov_count.in_R1])]:
for mn, exp in [(30, (0.09, 0.03, 0.04, 0.01)), (40, (0.19, 0.04, 0.06, 0.01))]:
d = d0[(d0.n >= mn) & d0.r.notna()]
a = d[d.genotype == "WT"]; b = d[d.genotype == "Fmr1-null"]
t = stats.ttest_ind(hc16.fisher(a.r), hc16.fisher(b.r))
u = stats.mannwhitneyu(a.r, b.r)
t1 = stats.ttest_ind(hc16.fisher(a.r_shuf1), hc16.fisher(b.r_shuf1))
t100 = stats.ttest_ind(a.r_shuf100, b.r_shuf100)
rows.append(dict(reading=label, min_passes=mn, n_WT=len(a), r_WT=a.r.mean(), sem_WT=a.r.sem(), n_KO=len(b),
r_KO=b.r.mean(), sem_KO=b.r.sem(), t=t.statistic, p_t=t.pvalue, U=u.statistic, p_U=u.pvalue,
d=hc16.cohen_d(hc16.fisher(a.r), hc16.fisher(b.r)),
shuf1_t=t1.statistic, shuf1_p=t1.pvalue, shuf100_t=t100.statistic, shuf100_p=t100.pvalue,
WT_in_tol=abs(a.r.mean() - exp[0]) <= exp[1], KO_in_tol=abs(b.r.mean() - exp[2]) <= exp[3]))
cv = pd.DataFrame(rows)
print(cv.round(4).to_string())
print("\nEach of the 100 permutations as a single-permutation test (t on Fisher z of the permuted r):")
for label, d0 in [("R2", cov), ("R2, count rule", cov_count)]:
for mn in [30, 40]:
d = d0[(d0.n >= mn) & d0.r.notna()]
a = np.vstack(d[d.genotype == "WT"].r_shuf_all); b = np.vstack(d[d.genotype == "Fmr1-null"].r_shuf_all)
tt = stats.ttest_ind(hc16.fisher(a), hc16.fisher(b), axis=0)
print(f" {label}, >= {mn}: p < 0.05 in {int((tt.pvalue < 0.05).sum())} of 100, p < 0.01 in {int((tt.pvalue < 0.01).sum())}; "
f"|t| median {np.median(np.abs(tt.statistic)):.2f}, max {np.abs(tt.statistic).max():.2f}; "
f"the first permutation (shuf1 above) t = {tt.statistic[0]:.2f}, ranked {int((np.abs(tt.statistic) >= abs(tt.statistic[0])).sum())} of 100 by |t|")
print("\nBy mouse (R2, >= 30 joint passes): pairs and mean r")
print(cov[cov.n >= 30].groupby(["genotype", "mouse"]).r.agg(["count", "mean"]).round(3).to_string())
print("\nBy mouse (R2, count rule, >= 30 joint passes): pairs and mean r")
print(cov_count[cov_count.n >= 30].groupby(["genotype", "mouse"]).r.agg(["count", "mean"]).round(3).to_string())
fig, ax = plt.subplots(figsize=(4, 3))
for g, col in [("WT", "k"), ("Fmr1-null", "r")]:
v = np.sort(cov[(cov.n >= 30) & (cov.genotype == g)].r.dropna())
ax.plot(v, np.arange(1, len(v) + 1) / len(v), color=col, label=g)
ax.set_xlabel("pair r over joint passes (>= 30)"); ax.set_ylabel("cumulative fraction"); ax.legend()
plt.show()
print("Plotted r range: %.3f to %.3f" % (cov[cov.n >= 30].r.min(), cov[cov.n >= 30].r.max()))
Joint passes per pair, literal rule: percentiles 0/10/50/90/100 = [39.0, 99.0, 177.0, 304.0, 346.0] ; pairs with >= 30: 332 and >= 40: 331 of 332
Joint passes per pair, count rule: percentiles 0/10/50/90/100 = [3.0, 9.0, 28.0, 73.0, 129.0] ; pairs with >= 30: 153 and >= 40: 116 of 328
Pairs with >= 3 joint passes: {'Fmr1-null': 210, 'WT': 122} ; joint passes range 39 to 346
reading min_passes n_WT r_WT sem_WT n_KO r_KO sem_KO t p_t U p_U d shuf1_t shuf1_p shuf100_t shuf100_p WT_in_tol KO_in_tol
0 R2 30 122 0.0608 0.0198 210 0.0210 0.0104 2.2890 0.0227 13752.0 0.2642 0.2606 1.2479 0.2130 1.2306 0.2193 True False
1 R2 40 122 0.0608 0.0198 209 0.0200 0.0104 2.3264 0.0206 13737.0 0.2397 0.2651 1.2082 0.2279 1.1428 0.2539 False False
2 R1 30 118 0.0611 0.0204 152 0.0170 0.0129 2.1393 0.0333 9754.0 0.2171 0.2625 0.7540 0.4515 0.8277 0.4086 True False
3 R1 40 118 0.0611 0.0204 151 0.0155 0.0129 2.1837 0.0299 9739.0 0.1902 0.2683 0.7049 0.4815 0.7151 0.4751 False False
4 R2, count rule 30 50 0.1646 0.0492 103 0.0074 0.0240 3.5426 0.0005 3212.0 0.0133 0.6106 3.6056 0.0004 0.6667 0.5060 False False
5 R2, count rule 40 30 0.2786 0.0662 86 0.0394 0.0256 4.4072 0.0000 1800.0 0.0013 0.9345 1.8708 0.0639 -0.0557 0.9557 False False
6 R1, count rule 30 48 0.1662 0.0509 63 0.0103 0.0331 2.8635 0.0050 1868.0 0.0343 0.5486 3.0207 0.0031 0.4445 0.6575 False False
7 R1, count rule 40 28 0.2894 0.0698 50 0.0480 0.0375 3.5108 0.0008 968.0 0.0053 0.8287 1.4744 0.1445 -0.3353 0.7383 False False
Each of the 100 permutations as a single-permutation test (t on Fisher z of the permuted r):
R2, >= 30: p < 0.05 in 4 of 100, p < 0.01 in 2; |t| median 0.67, max 3.22; the first permutation (shuf1 above) t = 1.25, ranked 19 of 100 by |t|
R2, >= 40: p < 0.05 in 3 of 100, p < 0.01 in 2; |t| median 0.67, max 3.18; the first permutation (shuf1 above) t = 1.21, ranked 20 of 100 by |t|
R2, count rule, >= 30: p < 0.05 in 4 of 100, p < 0.01 in 2; |t| median 0.68, max 3.61; the first permutation (shuf1 above) t = 3.61, ranked 1 of 100 by |t|
R2, count rule, >= 40: p < 0.05 in 5 of 100, p < 0.01 in 2; |t| median 0.71, max 2.86; the first permutation (shuf1 above) t = 1.87, ranked 8 of 100 by |t|
By mouse (R2, >= 30 joint passes): pairs and mean r
count mean
genotype mouse
Fmr1-null M18 22 0.036
M29 21 -0.001
M33 3 0.028
M34 102 0.036
M37 61 -0.005
M9 1 0.207
WT M1 6 0.067
M16 52 0.091
M28 7 -0.030
M32 45 0.049
M35 11 0.008
M36 1 0.205
By mouse (R2, count rule, >= 30 joint passes): pairs and mean r
count mean
genotype mouse
Fmr1-null M18 8 -0.013
M29 8 0.002
M33 1 0.040
M34 61 0.019
M37 24 -0.025
M9 1 0.261
WT M1 3 0.193
M16 18 0.252
M28 1 0.304
M32 20 0.150
M35 7 -0.065
M36 1 0.276
Plotted r range: -0.534 to 0.965
Result: partially reproduced, for both claims. fig4b-covariance-30-passes: in R2, 122 WT and 210 Fmr1-null pairs have at least 30 joint passes (expected 109 and 249). WT r = 0.061 +/- 0.020 and Fmr1-null 0.021 +/- 0.010; t on Fisher z = 2.29, p = 0.023 (Mann-Whitney p = 0.26). WT is within tolerance of 0.09 +/- 0.03, Fmr1-null is not within 0.04 +/- 0.01. After one permutation per pair the difference is gone (p = 0.21; mean of 100 permutations p = 0.22). R1 gives the same pattern (p = 0.033).
fig4b-covariance-40-passes, partially reproduced: under the literal pass rule the fewest joint passes of any pair is 39, so the 30 and 40 pass minimums select nothing and the sample hardly changes (122 WT, 209 Fmr1-null pairs): r = 0.061 against 0.020, t = 2.33, p = 0.021. The difference holds but is not larger, and both averages are outside tolerance (0.19 +/- 0.04, 0.06 +/- 0.01).
With the count rule pairs have 3 to 129 joint passes (median 28), and the minimums select: at 30, 50 WT and 103 Fmr1-null pairs, r = 0.165 against 0.007, t = 3.54, p = 0.0005; at 40, 30 and 86 pairs, r = 0.279 against 0.039, t = 4.41, p < 0.0001, a larger difference, as claimed. After 100 permutations neither differs (p = 0.51, 0.96), but one permutation per pair still differs at 30 passes (p = 0.0004). That draw is a chance one: taken one at a time, the 100 permutations give p < 0.05 in 4 of 100 at 30 passes (5 at 40), as a null should, and the first permutation, the one reported, has the largest |t| of the 100 (3.61; median 0.68). Under the literal rule 4 and 3 of 100 give p < 0.05. No average is within tolerance under the count rule (WT above, Fmr1-null below), so both outcomes would stay partially reproduced.
Kendall correlations¶
STAR Methods, "Functional Coupling of Cell Pairs": Kendall's correlation between the spike count series of each pair of simultaneously recorded cells, at 25, 40, 125, 250 ms, 1 s and 5 s. Tau-b, as in the CorrelateBPFSpikeTrainEpisodes code. The claims file names the readings: K1, every bin; K2, only bins in which the mouse is in an active pixel (rate at least the cell's mean, 32 x 32 map) of both cells; P1, the whole session; P2, the first 5 min; a pair counted once (R1 sessions) or per session (R2). Bins are placed on the spike clock from the first position sample; the position of a bin is the sample at its centre (hc16.Session.sample_time). Pairs are computed for every stationary session (R2) and the periods: the session, the first 5 min, the two halves and the second 5 min.
In [19]:
T_K = time.time()
def session_taus(t, cls, gate, extra_periods=True):
out = []
for sid, g in t[t.cls == cls].groupby("session"):
cl = list(g.cell)
if len(cl) < 2:
continue
s = S[sid]
D = s.sample_time[-1] - s.sample_time[0]
periods = {"session": (None, None), "first5": (0, 300)}
if extra_periods:
periods.update({"half1": (0, D / 2), "half2": (D / 2, D), "second5": (300, 600)})
out.append(hc16.pair_taus(s, cl, periods=periods, gate=gate))
d = pd.concat(out, ignore_index=True)
d["in_R1"] = d.session.isin(set(R1.session))
d["K"] = "K2" if gate else "K1"
return d
KT_place = pd.concat([session_taus(R2, "place", False), session_taus(R2, "place", True)], ignore_index=True)
KT_ns = pd.concat([session_taus(R2, "nonspatial", False), session_taus(R2, "nonspatial", True)], ignore_index=True)
print(f"Kendall computed in {time.time() - T_K:.0f} s")
for name, d in [("place", KT_place), ("non-spatial", KT_ns)]:
q = d[(d.period == "session") & (d.width == 0.025) & (d.K == "K1")]
print(f"{name} pairs (R2): {q.groupby('genotype').size().to_dict()}; in R1 sessions: {q[q.in_R1].groupby('genotype').size().to_dict()}")
print(" tau range over all readings %.3f to %.3f; undefined tau: %d of %d" % (d.tau.min(), d.tau.max(), d.tau.isna().sum(), len(d)))
print(" median bins entering K2 (session):", d[(d.K == "K2") & (d.period == "session")].groupby("width").nbins.median().to_dict())
def mw_table(d, expected=None):
rows = []
for (K, period, unit), _ in [((K, p, u), None) for K in ["K1", "K2"] for p in ["session", "first5"] for u in ["R1", "R2"]]:
q = d[(d.K == K) & (d.period == period)]
if unit == "R1":
q = q[q.in_R1]
for w in hc16.TIMESCALES:
qq = q[q.width == w].dropna(subset=["tau"])
a = qq[qq.genotype == "WT"].tau; b = qq[qq.genotype == "Fmr1-null"].tau
u = stats.mannwhitneyu(a, b)
row = dict(K=K, period=period, unit=unit, width=w, n_WT=len(a), WT=a.mean(), WT_sem=a.sem(), n_KO=len(b),
KO=b.mean(), KO_sem=b.sem(), U=u.statistic, p=u.pvalue, KO_lower=b.mean() < a.mean(),
d=hc16.cohen_d(a, b))
if expected and w in expected:
ew, sw, ek, sk = expected[w]
row["WT_in_tol"] = abs(a.mean() - ew) <= sw; row["KO_in_tol"] = abs(b.mean() - ek) <= sk
rows.append(row)
return pd.DataFrame(rows)
Kendall computed in 11 s
place pairs (R2): {'Fmr1-null': 210, 'WT': 122}; in R1 sessions: {'Fmr1-null': 152, 'WT': 118}
tau range over all readings -1.000 to 1.000; undefined tau: 845 of 19920
median bins entering K2 (session): {0.025: 3708.5, 0.04: 2325.0, 0.125: 741.5, 0.25: 366.5, 1.0: 93.5, 5.0: 18.0}
non-spatial pairs (R2): {'Fmr1-null': 146, 'WT': 380}; in R1 sessions: {'Fmr1-null': 115, 'WT': 205}
tau range over all readings -1.000 to 1.000; undefined tau: 2759 of 31560
median bins entering K2 (session): {0.025: 3960.5, 0.04: 2476.5, 0.125: 799.5, 0.25: 398.5, 1.0: 102.0, 5.0: 19.0}
fig4cd-kendall-place-cells¶
Kendall correlations between Fmr1-null place cell spike trains are weaker than between WT ones at the 25 ms, 40 ms and 125 ms timescales, but not significantly at 1 s and 5 s.
Criterion: Mann-Whitney between genotypes, pairs as units; reproduced when, under one reading, 25, 40 and 125 ms give p < 0.05 with Fmr1-null lower and 1 s and 5 s give p >= 0.05; partially reproduced when two of the three short timescales do. Pair counts expected 247 and 504.
In [20]:
exp_pc = {0.025: (0.017, 0.006, 0.009, 0.0012), 0.040: (0.02, 0.007, 0.01, 0.002), 0.125: (0.02, 0.008, 0.017, 0.003)}
kp = mw_table(KT_place, exp_pc)
print(kp.round(4).to_string())
def judge(kt, short, long_):
out = []
for key, g in kt.groupby(["K", "period", "unit"]):
g = g.set_index("width")
sig_short = [(g.loc[w, "p"] < 0.05) and g.loc[w, "KO_lower"] for w in short]
ns_long = [g.loc[w, "p"] >= 0.05 for w in long_]
out.append(dict(K=key[0], period=key[1], unit=key[2], short_sig=sum(sig_short), long_ns=sum(ns_long),
n_WT=int(g.n_WT.iloc[0]), n_KO=int(g.n_KO.iloc[0])))
return pd.DataFrame(out)
print("\nPer reading: short timescales with p < 0.05 and Fmr1-null lower; long timescales with p >= 0.05")
print(judge(kp, [0.025, 0.040, 0.125], [1.0, 5.0]).to_string())
print("\nBy mouse (K1, session, R2): pairs and mean tau at 25 ms")
q = KT_place[(KT_place.K == "K1") & (KT_place.period == "session") & (KT_place.width == 0.025)]
print(q.groupby(["genotype", "mouse"]).tau.agg(["count", "mean"]).round(4).to_string())
print("\nWithout M16 (52 of the 122 WT pairs): per reading, and K1 over the session")
kp16 = mw_table(KT_place[KT_place.mouse != "M16"], exp_pc)
print(judge(kp16, [0.025, 0.040, 0.125], [1.0, 5.0]).to_string())
print(kp16[(kp16.K == "K1") & (kp16.period == "session")].round(4).to_string())
print("Short timescales with p < 0.05 without M16, any reading:")
print(kp16[(kp16.width <= 0.125) & (kp16.p < 0.05)].round(4).to_string())
K period unit width n_WT WT WT_sem n_KO KO KO_sem U p KO_lower d WT_in_tol KO_in_tol
0 K1 session R1 0.025 118 0.0305 0.0122 152 0.0045 0.0024 9466.0 0.4344 True 0.2878 False False
1 K1 session R1 0.040 118 0.0334 0.0127 152 0.0057 0.0030 9229.0 0.6823 True 0.2909 False False
2 K1 session R1 0.125 118 0.0385 0.0143 152 0.0049 0.0054 9722.0 0.2364 True 0.2944 False False
3 K1 session R1 0.250 118 0.0411 0.0153 152 0.0048 0.0072 9806.0 0.1882 True 0.2819 NaN NaN
4 K1 session R1 1.000 118 0.0515 0.0176 152 0.0023 0.0118 10409.0 0.0236 True 0.2941 NaN NaN
5 K1 session R1 5.000 118 0.0801 0.0198 152 0.0009 0.0162 11046.0 0.0011 True 0.3833 NaN NaN
6 K1 session R2 0.025 122 0.0294 0.0118 210 0.0067 0.0022 12674.0 0.8723 True 0.2732 False False
7 K1 session R2 0.040 122 0.0321 0.0123 210 0.0084 0.0028 12447.0 0.6673 True 0.2693 False True
8 K1 session R2 0.125 122 0.0368 0.0138 210 0.0082 0.0046 13249.0 0.6030 True 0.2684 False False
9 K1 session R2 0.250 122 0.0394 0.0149 210 0.0079 0.0061 13495.0 0.4169 True 0.2584 NaN NaN
10 K1 session R2 1.000 122 0.0493 0.0171 210 0.0059 0.0099 14468.0 0.0493 True 0.2689 NaN NaN
11 K1 session R2 5.000 122 0.0770 0.0193 210 0.0032 0.0131 15535.0 0.0012 True 0.3712 NaN NaN
12 K1 first5 R1 0.025 118 0.0278 0.0124 152 0.0031 0.0021 8861.0 0.8671 True 0.2721 False False
13 K1 first5 R1 0.040 118 0.0296 0.0127 152 0.0024 0.0026 9037.0 0.9143 True 0.2879 False False
14 K1 first5 R1 0.125 118 0.0378 0.0147 152 0.0043 0.0051 9488.0 0.4144 True 0.2890 False False
15 K1 first5 R1 0.250 118 0.0410 0.0158 152 0.0054 0.0069 9495.0 0.4081 True 0.2735 NaN NaN
16 K1 first5 R1 1.000 118 0.0546 0.0178 152 0.0060 0.0110 10277.0 0.0398 True 0.2985 NaN NaN
17 K1 first5 R1 5.000 118 0.0836 0.0203 152 -0.0068 0.0162 11137.0 0.0007 True 0.4317 NaN NaN
18 K1 first5 R2 0.025 122 0.0267 0.0120 210 0.0052 0.0020 12197.0 0.4676 True 0.2586 False False
19 K1 first5 R2 0.040 122 0.0284 0.0123 210 0.0061 0.0025 12300.0 0.5457 True 0.2561 False False
20 K1 first5 R2 0.125 122 0.0356 0.0143 210 0.0077 0.0045 12924.0 0.8929 True 0.2570 False False
21 K1 first5 R2 0.250 122 0.0385 0.0153 210 0.0083 0.0060 13120.0 0.7136 True 0.2446 NaN NaN
22 K1 first5 R2 1.000 122 0.0519 0.0173 210 0.0093 0.0097 14482.0 0.0474 True 0.2645 NaN NaN
23 K1 first5 R2 5.000 122 0.0805 0.0199 210 -0.0039 0.0135 15664.0 0.0007 True 0.4131 NaN NaN
24 K2 session R1 0.025 118 0.0279 0.0137 149 -0.0000 0.0034 8987.0 0.7551 True 0.2698 False False
25 K2 session R1 0.040 118 0.0281 0.0144 149 0.0001 0.0044 8748.0 0.9459 True 0.2510 False False
26 K2 session R1 0.125 118 0.0249 0.0164 149 -0.0101 0.0076 9226.0 0.4881 True 0.2555 True False
27 K2 session R1 0.250 118 0.0130 0.0183 149 -0.0206 0.0106 9007.0 0.7309 True 0.2054 NaN NaN
28 K2 session R1 1.000 118 0.0144 0.0228 145 -0.0442 0.0154 9727.0 0.0562 True 0.2721 NaN NaN
29 K2 session R1 5.000 110 0.0396 0.0302 136 -0.0036 0.0270 7797.5 0.5678 True 0.1370 NaN NaN
30 K2 session R2 0.025 122 0.0268 0.0132 207 0.0042 0.0031 12163.0 0.5781 True 0.2363 False False
31 K2 session R2 0.040 122 0.0269 0.0140 207 0.0060 0.0039 11799.0 0.3207 True 0.2008 True False
32 K2 session R2 0.125 122 0.0236 0.0159 207 -0.0034 0.0064 12667.0 0.9622 True 0.2085 True False
33 K2 session R2 0.250 122 0.0120 0.0177 207 -0.0145 0.0086 12554.0 0.9307 True 0.1710 NaN NaN
34 K2 session R2 1.000 122 0.0132 0.0221 203 -0.0371 0.0123 13949.0 0.0563 True 0.2467 NaN NaN
35 K2 session R2 5.000 114 0.0449 0.0296 193 -0.0004 0.0205 11519.5 0.4906 True 0.1528 NaN NaN
36 K2 first5 R1 0.025 113 0.0224 0.0144 135 0.0021 0.0069 7767.0 0.8049 True 0.1701 True False
37 K2 first5 R1 0.040 113 0.0227 0.0153 135 -0.0016 0.0077 7979.0 0.5327 True 0.1901 True False
38 K2 first5 R1 0.125 113 0.0237 0.0192 132 0.0049 0.0138 7533.5 0.8921 True 0.1036 True False
39 K2 first5 R1 0.250 109 0.0084 0.0253 130 -0.0126 0.0194 7162.0 0.8857 True 0.0868 NaN NaN
40 K2 first5 R1 1.000 106 0.0489 0.0349 115 -0.0146 0.0287 6658.5 0.2358 True 0.1902 NaN NaN
41 K2 first5 R1 5.000 62 0.0399 0.0588 73 0.0042 0.0552 2361.5 0.6650 True 0.0764 NaN NaN
42 K2 first5 R2 0.025 117 0.0218 0.0140 191 0.0050 0.0053 10865.0 0.6847 True 0.1527 True False
43 K2 first5 R2 0.040 117 0.0220 0.0148 191 0.0054 0.0060 10889.0 0.7081 True 0.1397 True False
44 K2 first5 R2 0.125 117 0.0229 0.0187 188 0.0076 0.0106 10732.5 0.7235 True 0.0906 True False
45 K2 first5 R2 0.250 113 0.0091 0.0246 186 -0.0113 0.0146 10490.0 0.9796 True 0.0910 NaN NaN
46 K2 first5 R2 1.000 110 0.0426 0.0338 170 -0.0084 0.0218 9933.5 0.3783 True 0.1624 NaN NaN
47 K2 first5 R2 5.000 65 0.0511 0.0567 122 0.0219 0.0396 4077.0 0.7516 True 0.0657 NaN NaN
Per reading: short timescales with p < 0.05 and Fmr1-null lower; long timescales with p >= 0.05
K period unit short_sig long_ns n_WT n_KO
0 K1 first5 R1 0 0 118 152
1 K1 first5 R2 0 0 122 210
2 K1 session R1 0 0 118 152
3 K1 session R2 0 0 122 210
4 K2 first5 R1 0 2 113 135
5 K2 first5 R2 0 2 117 191
6 K2 session R1 0 2 118 149
7 K2 session R2 0 2 122 207
By mouse (K1, session, R2): pairs and mean tau at 25 ms
count mean
genotype mouse
Fmr1-null M18 22 0.0182
M29 21 0.0080
M33 3 0.0074
M34 102 0.0083
M37 61 -0.0011
M9 1 0.0349
WT M1 6 0.0220
M16 52 0.0522
M28 7 0.0020
M32 45 0.0134
M35 11 0.0100
M36 1 0.0062
Without M16 (52 of the 122 WT pairs): per reading, and K1 over the session
K period unit short_sig long_ns n_WT n_KO
0 K1 first5 R1 1 0 67 152
1 K1 first5 R2 0 0 70 210
2 K1 session R1 0 0 67 152
3 K1 session R2 0 0 70 210
4 K2 first5 R1 0 2 63 135
5 K2 first5 R2 0 2 66 191
6 K2 session R1 0 2 67 149
7 K2 session R2 0 2 70 207
K period unit width n_WT WT WT_sem n_KO KO KO_sem U p KO_lower d WT_in_tol KO_in_tol
0 K1 session R1 0.025 67 0.0131 0.0053 152 0.0045 0.0024 5408.0 0.4653 True 0.2524 True False
1 K1 session R1 0.040 67 0.0171 0.0068 152 0.0057 0.0030 5383.0 0.5014 True 0.2608 True False
2 K1 session R1 0.125 67 0.0259 0.0107 152 0.0049 0.0054 5824.0 0.0905 True 0.2860 True False
3 K1 session R1 0.250 67 0.0326 0.0134 152 0.0048 0.0072 5906.0 0.0597 True 0.2905 NaN NaN
4 K1 session R1 1.000 67 0.0525 0.0192 152 0.0023 0.0118 6191.0 0.0110 True 0.3365 NaN NaN
5 K1 session R1 5.000 67 0.0775 0.0241 152 0.0009 0.0162 6327.0 0.0043 True 0.3855 NaN NaN
6 K1 session R2 0.025 70 0.0124 0.0051 210 0.0067 0.0022 7334.0 0.9789 True 0.1630 True False
7 K1 session R2 0.040 70 0.0161 0.0065 210 0.0084 0.0028 7353.0 0.9966 True 0.1751 True True
8 K1 session R2 0.125 70 0.0242 0.0103 210 0.0082 0.0046 8001.0 0.2675 True 0.2214 True False
9 K1 session R2 0.250 70 0.0303 0.0129 210 0.0079 0.0061 8180.0 0.1574 True 0.2389 NaN NaN
10 K1 session R2 1.000 70 0.0482 0.0185 210 0.0059 0.0099 8662.0 0.0254 True 0.2900 NaN NaN
11 K1 session R2 5.000 70 0.0728 0.0233 210 0.0032 0.0131 8967.0 0.0059 True 0.3637 NaN NaN
Short timescales with p < 0.05 without M16, any reading:
K period unit width n_WT WT WT_sem n_KO KO KO_sem U p KO_lower d WT_in_tol KO_in_tol
14 K1 first5 R1 0.125 67 0.0269 0.0111 152 0.0043 0.0051 5957.0 0.0454 True 0.3109 True False
Result: not reproduced. Under no reading do the 25, 40 and 125 ms comparisons reach p < 0.05 (smallest p over the eight readings 0.24). With K1 over the session (R2, 122 WT and 210 Fmr1-null pairs, expected 247 and 504) mean tau at 25 ms is 0.029 against 0.0067, lower in Fmr1-null, but Mann-Whitney p = 0.87, 0.67 and 0.60 at 25, 40 and 125 ms, while 1 s (p = 0.049) and 5 s (p = 0.0012) do differ, the reverse of the claimed pattern. No mean tau at 25 ms is within tolerance under the K1 readings. Without M16, which gives 52 of the 122 WT pairs, the pattern is the same (K1 over the session, R2: p = 0.98, 1.00 and 0.27 at 25, 40 and 125 ms; 0.025 and 0.0059 at 1 s and 5 s), so it is not one mouse's; one of the eight readings without M16 has one short timescale below 0.05 (K1, first 5 min, R1, 125 ms, p = 0.045).
fig4d-network-state-recurrence¶
Network states recur less in Fmr1-null ensembles; the vector of pairwise Kendall correlations does not differ between the two halves of WT recordings but does in Fmr1-null recordings.
Criterion: for each pair, tau in the earlier and later interval; a paired t test per genotype and timescale, and PCo, the Pearson correlation of the two vectors (all pairs of a genotype pooled). Readings: H1 the two halves, H2 the first two 5 min intervals, each with K1 and K2. Reproduced when, under one reading, no WT test has p < 0.05, at least five of the six Fmr1-null tests do, and WT PCo exceeds Fmr1-null PCo at most timescales. Pairs per session (R2), as the paper's 247 and 504 pairs; R1 shown too.
In [21]:
rows = []
for K in ["K1", "K2"]:
for H, (p1, p2) in [("H1", ("half1", "half2")), ("H2", ("first5", "second5"))]:
for unit in ["R2", "R1"]:
for g in ["WT", "Fmr1-null"]:
for w in hc16.TIMESCALES:
q = KT_place[(KT_place.K == K) & (KT_place.genotype == g) & (KT_place.width == w)]
if unit == "R1":
q = q[q.in_R1]
a = q[q.period == p1].set_index(["session", "a", "b"]).tau
b = q[q.period == p2].set_index(["session", "a", "b"]).tau
j = pd.concat([a, b], axis=1, keys=["x", "y"]).dropna()
t = stats.ttest_rel(j.x, j.y)
rows.append(dict(K=K, H=H, unit=unit, genotype=g, width=w, n=len(j), t=t.statistic, p=t.pvalue,
PCo=stats.pearsonr(j.x, j.y).statistic))
rec = pd.DataFrame(rows)
print(rec.round(4).to_string())
summ = []
for key, g in rec.groupby(["K", "H", "unit"]):
wt = g[g.genotype == "WT"].set_index("width"); ko = g[g.genotype == "Fmr1-null"].set_index("width")
summ.append(dict(K=key[0], H=key[1], unit=key[2], WT_sig=int((wt.p < 0.05).sum()), KO_sig=int((ko.p < 0.05).sum()),
PCo_WT_gt_KO=int((wt.PCo > ko.PCo).sum())))
print("\nPer reading (of 6 timescales):")
print(pd.DataFrame(summ).to_string())
K H unit genotype width n t p PCo
0 K1 H1 R2 WT 0.025 122 0.1997 0.8420 0.9908
1 K1 H1 R2 WT 0.040 122 0.4651 0.6427 0.9869
2 K1 H1 R2 WT 0.125 122 0.7349 0.4638 0.9676
3 K1 H1 R2 WT 0.250 122 0.6790 0.4984 0.9486
4 K1 H1 R2 WT 1.000 122 1.0543 0.2938 0.8950
5 K1 H1 R2 WT 5.000 122 0.3347 0.7384 0.8028
6 K1 H1 R2 Fmr1-null 0.025 210 -1.2727 0.2045 0.8265
7 K1 H1 R2 Fmr1-null 0.040 210 -1.1861 0.2369 0.8112
8 K1 H1 R2 Fmr1-null 0.125 210 -0.5764 0.5650 0.8397
9 K1 H1 R2 Fmr1-null 0.250 210 -0.3666 0.7143 0.8308
10 K1 H1 R2 Fmr1-null 1.000 210 -0.1663 0.8681 0.8333
11 K1 H1 R2 Fmr1-null 5.000 210 0.0927 0.9263 0.7762
12 K1 H1 R1 WT 0.025 118 0.1721 0.8637 0.9908
13 K1 H1 R1 WT 0.040 118 0.4068 0.6849 0.9869
14 K1 H1 R1 WT 0.125 118 0.6963 0.4877 0.9676
15 K1 H1 R1 WT 0.250 118 0.6803 0.4977 0.9485
16 K1 H1 R1 WT 1.000 118 1.1127 0.2681 0.8947
17 K1 H1 R1 WT 5.000 118 0.4084 0.6837 0.8015
18 K1 H1 R1 Fmr1-null 0.025 152 -0.8956 0.3719 0.7534
19 K1 H1 R1 Fmr1-null 0.040 152 -0.9714 0.3329 0.7441
20 K1 H1 R1 Fmr1-null 0.125 152 -0.3855 0.7004 0.7998
21 K1 H1 R1 Fmr1-null 0.250 152 -0.2747 0.7839 0.7961
22 K1 H1 R1 Fmr1-null 1.000 152 -0.5404 0.5897 0.8133
23 K1 H1 R1 Fmr1-null 5.000 152 -0.4957 0.6208 0.7599
24 K1 H2 R2 WT 0.025 122 -2.1317 0.0351 0.9838
25 K1 H2 R2 WT 0.040 122 -2.4362 0.0163 0.9792
26 K1 H2 R2 WT 0.125 122 -1.4795 0.1416 0.9485
27 K1 H2 R2 WT 0.250 122 -1.2162 0.2263 0.9197
28 K1 H2 R2 WT 1.000 122 -0.8149 0.4167 0.8416
29 K1 H2 R2 WT 5.000 122 0.4317 0.6667 0.7206
30 K1 H2 R2 Fmr1-null 0.025 207 0.1115 0.9114 0.7922
31 K1 H2 R2 Fmr1-null 0.040 207 -0.6728 0.5018 0.8191
32 K1 H2 R2 Fmr1-null 0.125 207 0.6888 0.4917 0.8225
33 K1 H2 R2 Fmr1-null 0.250 207 0.3855 0.7003 0.8159
34 K1 H2 R2 Fmr1-null 1.000 207 1.3918 0.1655 0.7876
35 K1 H2 R2 Fmr1-null 5.000 207 -0.5097 0.6108 0.6956
36 K1 H2 R1 WT 0.025 118 -2.0365 0.0440 0.9838
37 K1 H2 R1 WT 0.040 118 -2.3539 0.0202 0.9793
38 K1 H2 R1 WT 0.125 118 -1.2222 0.2241 0.9496
39 K1 H2 R1 WT 0.250 118 -0.9371 0.3506 0.9217
40 K1 H2 R1 WT 1.000 118 -0.7570 0.4506 0.8426
41 K1 H2 R1 WT 5.000 118 0.4916 0.6239 0.7222
42 K1 H2 R1 Fmr1-null 0.025 149 0.3171 0.7517 0.6905
43 K1 H2 R1 Fmr1-null 0.040 149 -1.1348 0.2583 0.7538
44 K1 H2 R1 Fmr1-null 0.125 149 0.6353 0.5262 0.7927
45 K1 H2 R1 Fmr1-null 0.250 149 0.7069 0.4808 0.7999
46 K1 H2 R1 Fmr1-null 1.000 149 1.7489 0.0824 0.7538
47 K1 H2 R1 Fmr1-null 5.000 149 -0.0275 0.9781 0.6965
48 K2 H1 R2 WT 0.025 121 0.5386 0.5911 0.9420
49 K2 H1 R2 WT 0.040 121 -0.7763 0.4391 0.8289
50 K2 H1 R2 WT 0.125 121 0.9055 0.3670 0.8404
51 K2 H1 R2 WT 0.250 120 -0.6057 0.5459 0.7941
52 K2 H1 R2 WT 1.000 118 0.9693 0.3344 0.6861
53 K2 H1 R2 WT 5.000 96 0.6923 0.4904 0.2138
54 K2 H1 R2 Fmr1-null 0.025 204 0.1499 0.8810 0.4684
55 K2 H1 R2 Fmr1-null 0.040 204 0.6202 0.5358 0.4086
56 K2 H1 R2 Fmr1-null 0.125 203 0.4601 0.6459 0.4350
57 K2 H1 R2 Fmr1-null 0.250 199 0.2419 0.8091 0.3681
58 K2 H1 R2 Fmr1-null 1.000 193 1.1866 0.2369 0.3106
59 K2 H1 R2 Fmr1-null 5.000 162 1.2260 0.2220 0.1904
60 K2 H1 R1 WT 0.025 117 0.4116 0.6814 0.9424
61 K2 H1 R1 WT 0.040 117 -0.8241 0.4116 0.8288
62 K2 H1 R1 WT 0.125 117 0.8226 0.4124 0.8427
63 K2 H1 R1 WT 0.250 116 -0.6741 0.5016 0.7946
64 K2 H1 R1 WT 1.000 114 1.0829 0.2811 0.6871
65 K2 H1 R1 WT 5.000 92 0.4635 0.6441 0.2056
66 K2 H1 R1 Fmr1-null 0.025 146 0.6455 0.5197 0.3255
67 K2 H1 R1 Fmr1-null 0.040 146 0.6968 0.4870 0.2518
68 K2 H1 R1 Fmr1-null 0.125 145 0.5401 0.5900 0.3114
69 K2 H1 R1 Fmr1-null 0.250 141 0.2439 0.8077 0.2756
70 K2 H1 R1 Fmr1-null 1.000 136 0.8665 0.3877 0.2413
71 K2 H1 R1 Fmr1-null 5.000 108 0.8671 0.3878 0.1662
72 K2 H2 R2 WT 0.025 116 -1.1203 0.2649 0.9099
73 K2 H2 R2 WT 0.040 115 -0.6088 0.5439 0.8761
74 K2 H2 R2 WT 0.125 116 -0.5366 0.5926 0.7743
75 K2 H2 R2 WT 0.250 113 0.3256 0.7453 0.6799
76 K2 H2 R2 WT 1.000 109 -0.2097 0.8343 0.2657
77 K2 H2 R2 WT 5.000 53 0.4950 0.6227 0.4631
78 K2 H2 R2 Fmr1-null 0.025 175 -0.1732 0.8627 0.3467
79 K2 H2 R2 Fmr1-null 0.040 175 -1.1205 0.2640 0.6682
80 K2 H2 R2 Fmr1-null 0.125 174 -0.0854 0.9320 0.6358
81 K2 H2 R2 Fmr1-null 0.250 168 -0.8567 0.3928 0.5921
82 K2 H2 R2 Fmr1-null 1.000 160 1.8090 0.0723 0.4991
83 K2 H2 R2 Fmr1-null 5.000 99 -0.1583 0.8745 0.0175
84 K2 H2 R1 WT 0.025 112 -0.9093 0.3652 0.9116
85 K2 H2 R1 WT 0.040 111 -0.4024 0.6882 0.8776
86 K2 H2 R1 WT 0.125 112 -0.4566 0.6488 0.7742
87 K2 H2 R1 WT 0.250 109 0.3624 0.7177 0.6837
88 K2 H2 R1 WT 1.000 105 0.0323 0.9743 0.2794
89 K2 H2 R1 WT 5.000 50 0.2232 0.8243 0.4685
90 K2 H2 R1 Fmr1-null 0.025 119 0.0213 0.9831 0.2897
91 K2 H2 R1 Fmr1-null 0.040 119 -1.4834 0.1406 0.6884
92 K2 H2 R1 Fmr1-null 0.125 118 0.4845 0.6289 0.6750
93 K2 H2 R1 Fmr1-null 0.250 112 -0.1987 0.8429 0.6221
94 K2 H2 R1 Fmr1-null 1.000 105 2.1571 0.0333 0.5193
95 K2 H2 R1 Fmr1-null 5.000 53 -0.2237 0.8239 -0.0606
Per reading (of 6 timescales):
K H unit WT_sig KO_sig PCo_WT_gt_KO
0 K1 H1 R1 0 0 6
1 K1 H1 R2 0 0 6
2 K1 H2 R1 2 0 6
3 K1 H2 R2 2 0 6
4 K2 H1 R1 0 0 6
5 K2 H1 R2 0 0 6
6 K2 H2 R1 0 1 5
7 K2 H2 R2 0 0 5
Result: not reproduced. No reading has at least five significant Fmr1-null tests: under H1 none of the six Fmr1-null paired tests reaches p < 0.05 with K1 or K2, and under H2 at most one. WT halves do not differ under H1; under H2 with K1 the WT first two 5 min intervals differ at 25 and 40 ms (p = 0.035, 0.016). PCo is higher in WT than in Fmr1-null at all six timescales in six of the eight readings (K1, H1, R2, 25 ms: 0.991 against 0.827) and at five in the other two (K2, H2, R1 and R2), where Fmr1-null is the higher at 1 s (0.519 against 0.279, 0.499 against 0.266).
kendall-nonspatial-pyramidal¶
Kendall correlations between Fmr1-null non-spatial pyramidal cells are weaker than between WT ones at 25 ms, 40 ms, 125 ms and 1 s, but not significantly at 5 s.
Criterion: as fig4cd for non-spatial pyramidal pairs; reproduced when, under one reading, 25, 40, 125 ms and 1 s give p < 0.05 with Fmr1-null lower and 5 s gives p >= 0.05; partially reproduced when three of the four do. Expected pairs 252 and 203.
In [22]:
exp_ns = {0.025: (0.03, 0.006, 0.007, 0.002), 0.125: (0.04, 0.007, 0.01, 0.003)}
kn = mw_table(KT_ns, exp_ns)
print(kn.round(4).to_string())
print("\nPer reading:")
print(judge(kn, [0.025, 0.040, 0.125, 1.0], [5.0]).to_string())
print("\nBy mouse (K1, session, R2): pairs and mean tau at 25 ms")
q = KT_ns[(KT_ns.K == "K1") & (KT_ns.period == "session") & (KT_ns.width == 0.025)]
print(q.groupby(["genotype", "mouse"]).tau.agg(["count", "mean"]).round(4).to_string())
K period unit width n_WT WT WT_sem n_KO KO KO_sem U p KO_lower d WT_in_tol KO_in_tol
0 K1 session R1 0.025 205 0.0083 0.0024 115 0.0110 0.0039 9142.0 0.0009 False -0.0737 False False
1 K1 session R1 0.040 205 0.0099 0.0027 115 0.0128 0.0040 9099.0 0.0007 False -0.0740 NaN NaN
2 K1 session R1 0.125 205 0.0139 0.0035 115 0.0170 0.0048 9430.0 0.0030 False -0.0622 False False
3 K1 session R1 0.250 205 0.0164 0.0041 115 0.0203 0.0053 9649.0 0.0071 False -0.0663 NaN NaN
4 K1 session R1 1.000 205 0.0215 0.0053 115 0.0263 0.0073 10269.0 0.0559 False -0.0635 NaN NaN
5 K1 session R1 5.000 205 0.0296 0.0068 115 0.0281 0.0109 11469.0 0.6888 True 0.0146 NaN NaN
6 K1 session R2 0.025 380 0.0089 0.0016 146 0.0099 0.0032 23626.0 0.0084 False -0.0306 False False
7 K1 session R2 0.040 380 0.0109 0.0018 146 0.0118 0.0033 22994.0 0.0024 False -0.0253 NaN NaN
8 K1 session R2 0.125 380 0.0152 0.0025 146 0.0155 0.0042 23396.0 0.0054 False -0.0063 False False
9 K1 session R2 0.250 380 0.0185 0.0030 146 0.0179 0.0050 24099.0 0.0197 True 0.0104 NaN NaN
10 K1 session R2 1.000 380 0.0243 0.0041 146 0.0207 0.0086 25288.0 0.1163 True 0.0407 NaN NaN
11 K1 session R2 5.000 380 0.0276 0.0057 146 0.0219 0.0136 26797.0 0.5460 True 0.0441 NaN NaN
12 K1 first5 R1 0.025 205 0.0097 0.0027 115 0.0095 0.0043 10986.0 0.3131 True 0.0053 False False
13 K1 first5 R1 0.040 205 0.0122 0.0032 115 0.0097 0.0045 11887.0 0.9008 True 0.0555 NaN NaN
14 K1 first5 R1 0.125 205 0.0159 0.0043 115 0.0153 0.0061 11710.5 0.9233 True 0.0098 False False
15 K1 first5 R1 0.250 205 0.0186 0.0051 115 0.0167 0.0074 12087.5 0.7061 True 0.0258 NaN NaN
16 K1 first5 R1 1.000 205 0.0249 0.0070 115 0.0125 0.0104 12669.0 0.2673 True 0.1186 NaN NaN
17 K1 first5 R1 5.000 205 0.0411 0.0101 115 -0.0016 0.0168 14097.0 0.0036 True 0.2698 NaN NaN
18 K1 first5 R2 0.025 380 0.0082 0.0017 146 0.0092 0.0035 26370.0 0.3803 False -0.0294 False False
19 K1 first5 R2 0.040 380 0.0106 0.0020 146 0.0102 0.0037 28139.0 0.7985 True 0.0117 NaN NaN
20 K1 first5 R2 0.125 380 0.0148 0.0029 146 0.0150 0.0052 27580.5 0.9189 False -0.0045 False False
21 K1 first5 R2 0.250 380 0.0182 0.0034 146 0.0165 0.0062 28219.5 0.7589 True 0.0233 NaN NaN
22 K1 first5 R2 1.000 380 0.0260 0.0049 146 0.0139 0.0100 29791.0 0.1890 True 0.1190 NaN NaN
23 K1 first5 R2 5.000 380 0.0362 0.0075 146 -0.0014 0.0155 32162.0 0.0046 True 0.2365 NaN NaN
24 K2 session R1 0.025 197 0.0274 0.0072 115 0.0236 0.0070 8998.0 0.0024 True 0.0403 True False
25 K2 session R1 0.040 197 0.0314 0.0078 115 0.0321 0.0086 9212.0 0.0059 False -0.0070 NaN NaN
26 K2 session R1 0.125 193 0.0419 0.0109 115 0.0445 0.0113 9603.0 0.0481 False -0.0188 True False
27 K2 session R1 0.250 193 0.0459 0.0132 115 0.0497 0.0130 9495.0 0.0341 False -0.0225 NaN NaN
28 K2 session R1 1.000 181 0.0570 0.0186 112 0.0629 0.0154 9120.5 0.1498 False -0.0265 NaN NaN
29 K2 session R1 5.000 155 0.0293 0.0254 107 0.1200 0.0239 6794.0 0.0130 False -0.3127 NaN NaN
30 K2 session R2 0.025 369 0.0251 0.0045 146 0.0214 0.0057 22680.0 0.0052 True 0.0449 True False
31 K2 session R2 0.040 368 0.0297 0.0051 146 0.0299 0.0071 22340.0 0.0029 False -0.0021 NaN NaN
32 K2 session R2 0.125 364 0.0392 0.0071 146 0.0411 0.0095 23155.5 0.0232 False -0.0143 True False
33 K2 session R2 0.250 364 0.0443 0.0086 146 0.0431 0.0113 24037.5 0.0921 True 0.0074 NaN NaN
34 K2 session R2 1.000 345 0.0638 0.0124 143 0.0505 0.0143 23616.5 0.4588 True 0.0618 NaN NaN
35 K2 session R2 5.000 310 0.0383 0.0172 138 0.1068 0.0230 18214.5 0.0121 False -0.2334 NaN NaN
36 K2 first5 R1 0.025 159 0.0316 0.0088 98 0.0103 0.0048 7397.0 0.4966 True 0.2311 True False
37 K2 first5 R1 0.040 156 0.0362 0.0093 99 0.0124 0.0054 7692.0 0.9590 True 0.2465 NaN NaN
38 K2 first5 R1 0.125 151 0.0432 0.0150 95 0.0272 0.0121 6819.0 0.5159 True 0.0984 True False
39 K2 first5 R1 0.250 155 0.0620 0.0188 95 0.0396 0.0182 7369.5 0.9907 True 0.1045 NaN NaN
40 K2 first5 R1 1.000 136 0.0746 0.0277 86 0.0457 0.0282 5907.0 0.9001 True 0.0959 NaN NaN
41 K2 first5 R1 5.000 76 0.0866 0.0553 68 0.0925 0.0490 2610.5 0.9171 False -0.0130 NaN NaN
42 K2 first5 R2 0.025 304 0.0253 0.0056 127 0.0155 0.0055 18140.0 0.3237 True 0.1114 True False
43 K2 first5 R2 0.040 302 0.0288 0.0060 128 0.0244 0.0084 18953.0 0.7506 True 0.0440 NaN NaN
44 K2 first5 R2 0.125 290 0.0382 0.0094 123 0.0299 0.0112 17248.0 0.5970 True 0.0554 True False
45 K2 first5 R2 0.250 285 0.0492 0.0123 123 0.0457 0.0176 17063.0 0.6712 True 0.0171 NaN NaN
46 K2 first5 R2 1.000 240 0.0754 0.0189 115 0.0662 0.0257 13701.0 0.9133 True 0.0320 NaN NaN
47 K2 first5 R2 5.000 166 0.0198 0.0363 89 0.0629 0.0440 7053.5 0.5530 False -0.0957 NaN NaN
Per reading:
K period unit short_sig long_ns n_WT n_KO
0 K1 first5 R1 0 0 205 115
1 K1 first5 R2 0 0 380 146
2 K1 session R1 0 1 205 115
3 K1 session R2 0 1 380 146
4 K2 first5 R1 0 1 159 98
5 K2 first5 R2 0 1 304 127
6 K2 session R1 1 0 197 115
7 K2 session R2 1 0 369 146
By mouse (K1, session, R2): pairs and mean tau at 25 ms
count mean
genotype mouse
Fmr1-null M10 1 0.0003
M18 35 0.0037
M29 1 0.0061
M33 21 0.0249
M34 32 0.0071
M37 53 0.0105
M9 3 -0.0000
WT M11 1 -0.0023
M16 297 0.0089
M28 10 0.0110
M32 35 0.0142
M35 22 0.0051
M36 15 0.0018
Result: not reproduced. There are 380 WT and 146 Fmr1-null non-spatial pairs in R2 (expected 252 and 203). With K1 over the session the distributions differ at 25, 40 and 125 ms (p = 0.0084, 0.0024, 0.0054) but Fmr1-null mean tau is not lower (0.0099 against 0.0089 at 25 ms), and 1 s does not differ (p = 0.12). No reading gives Fmr1-null lower at p < 0.05 on more than one of the four timescales (K2 over the session, 25 ms, p = 0.0024 in R1).
Refractory filter and 32 x 32 grid for the pass and Kendall claims¶
The two variants of the cell table (refractory filter, 32 x 32 grid; see the cell table section) applied to the claims that use passes or Kendall correlations, R2 only. With the refractory filter a pass or a pair enters only when its cells pass the filter. On the 32 x 32 grid the place cells are those classified on the grid, and exp is computed from their 32 x 32 rate maps. The covariance shuffle uses its own generator. fig4d-network-state-recurrence, kendall-quiet-states and the rotating arena claims are not run under these variants.
In [23]:
T_V = time.time()
keep_iso = set(zip(R2_iso.session, R2_iso.cell))
def kept(d, *cols):
k = np.ones(len(d), bool)
for col in cols:
k &= np.array([(s, c) in keep_iso for s, c in zip(d.session, d[col])])
return d[k]
def zframe(P):
return pd.DataFrame([dict(session=k[0], cell=k[1], genotype=S[k[0]].genotype, z=zz)
for k, v in P.items() for zz in v["z"][v["is_pass"]]])
P_code = pass_table(R2_code, index_fn=hc16.pixel_index_code)
P_code_count = pass_table(R2_code, rule="count", index_fn=hc16.pixel_index_code)
variants = {"literal": {"as above": pz, "refractory filter": kept(pz, "cell"), "32 x 32 grid": zframe(P_code)},
"count": {"as above": pz_count, "refractory filter": kept(pz_count, "cell"), "32 x 32 grid": zframe(P_code_count)}}
print("R2 place cells with passes:", {f"{rule}, {name}": zz.groupby("genotype")[["session", "cell"]].apply(lambda d: len(d.drop_duplicates())).to_dict()
for rule, vv in variants.items() for name, zz in vv.items()})
print("\noverdispersion-passes (R2), variance of z over passes:")
print(pd.DataFrame([od_passes(zz, f"{rule} rule, {name}") for rule, vv in variants.items() for name, zz in vv.items()]).round(4).to_string())
print("\noverdispersion-cells (R2), mean of the per-cell variance of z, one-tailed WT > Fmr1-null:")
rows = []
for rule, vv in variants.items():
for name, zz in vv.items():
pv = zz.groupby(["session", "cell", "genotype"]).z.var().reset_index()
r = hc16.two_group(pv[pv.genotype == "WT"].z, pv[pv.genotype == "Fmr1-null"].z, alternative="greater")
rows.append(dict(reading=f"{rule} rule, {name}", n_WT=r["n_a"], WT=r["mean_a"], WT_sem=r["sem_a"], n_KO=r["n_b"], KO=r["mean_b"],
KO_sem=r["sem_b"], p_t_one=r["p_t"], p_U_one=r["p_U"],
WT_in_tol=abs(r["mean_a"] - 5.98) <= 0.44, KO_in_tol=abs(r["mean_b"] - 5.02) <= 0.3))
print(pd.DataFrame(rows).round(4).to_string())
print("\nfig4b-covariance-30-passes and -40-passes (R2), pair r over joint passes:")
rng_var = np.random.default_rng(406)
covs = {"literal": {"as above": cov, "refractory filter": kept(cov, "a", "b"), "32 x 32 grid": pair_cov(P_code, rng=rng_var)},
"count": {"as above": cov_count, "refractory filter": kept(cov_count, "a", "b"), "32 x 32 grid": pair_cov(P_code_count, rng=rng_var)}}
rows = []
for rule, vv in covs.items():
for name, d0 in vv.items():
for mn, exp in [(30, (0.09, 0.03, 0.04, 0.01)), (40, (0.19, 0.04, 0.06, 0.01))]:
d = d0[(d0.n >= mn) & d0.r.notna()]
a = d[d.genotype == "WT"]; b = d[d.genotype == "Fmr1-null"]
t = stats.ttest_ind(hc16.fisher(a.r), hc16.fisher(b.r))
rows.append(dict(reading=f"{rule} rule, {name}", min_passes=mn, joint_min=int(d0.n.min()), n_WT=len(a), r_WT=a.r.mean(),
n_KO=len(b), r_KO=b.r.mean(), t=t.statistic, p_t=t.pvalue,
shuf1_p=stats.ttest_ind(hc16.fisher(a.r_shuf1), hc16.fisher(b.r_shuf1)).pvalue,
shuf100_p=stats.ttest_ind(a.r_shuf100, b.r_shuf100).pvalue,
WT_in_tol=abs(a.r.mean() - exp[0]) <= exp[1], KO_in_tol=abs(b.r.mean() - exp[2]) <= exp[3]))
print(pd.DataFrame(rows).round(4).to_string())
KT_place_code = pd.concat([session_taus(R2_code, "place", False), session_taus(R2_code, "place", True)], ignore_index=True)
KT_ns_code = pd.concat([session_taus(R2_code, "nonspatial", False), session_taus(R2_code, "nonspatial", True)], ignore_index=True)
for claim, kts, short, long_ in [("fig4cd-kendall-place-cells", {"refractory filter": kept(KT_place, "a", "b"), "32 x 32 grid": KT_place_code},
[0.025, 0.040, 0.125], [1.0, 5.0]),
("kendall-nonspatial-pyramidal", {"refractory filter": kept(KT_ns, "a", "b"), "32 x 32 grid": KT_ns_code},
[0.025, 0.040, 0.125, 1.0], [5.0])]:
for name, kt in kts.items():
m = mw_table(kt)
print(f"\n{claim}, {name}: per reading, timescales {short} with p < 0.05 and Fmr1-null lower, and {long_} with p >= 0.05")
print(judge(m, short, long_).to_string())
print(" K1, session, R2:", "; ".join(f"{w * 1000:g} ms: WT {r.WT:.4f}, KO {r.KO:.4f}, p = {r.p:.3g}"
for w, r in m[(m.K == "K1") & (m.period == "session") & (m.unit == "R2")].set_index("width").iterrows()))
print(f"\nVariants computed in {time.time() - T_V:.0f} s")
R2 place cells with passes: {'literal, as above': {'Fmr1-null': 70, 'WT': 55}, 'literal, refractory filter': {'Fmr1-null': 70, 'WT': 54}, 'literal, 32 x 32 grid': {'Fmr1-null': 81, 'WT': 65}, 'count, as above': {'Fmr1-null': 70, 'WT': 55}, 'count, refractory filter': {'Fmr1-null': 70, 'WT': 54}, 'count, 32 x 32 grid': {'Fmr1-null': 81, 'WT': 65}}
overdispersion-passes (R2), variance of z over passes:
reading n_WT var_WT n_KO var_KO F p_F_one_sided p_F_two_sided levene_F levene_p WT_in_tol KO_in_tol passes_in_tol
0 literal rule, as above 12802 3.6903 18521 3.5493 1.0397 0.0082 0.0164 0.0000 0.9977 False False False
1 literal rule, refractory filter 12652 3.7075 18521 3.5493 1.0446 0.0037 0.0074 0.0298 0.8630 False False False
2 literal rule, 32 x 32 grid 15777 3.7305 22036 3.7297 1.0002 0.4940 0.9880 4.7290 0.0297 False False False
3 count rule, as above 5442 5.2789 8250 4.6635 1.1319 0.0000 0.0000 11.6845 0.0006 True True False
4 count rule, refractory filter 5386 5.3018 8250 4.6635 1.1369 0.0000 0.0000 12.5852 0.0004 True True False
5 count rule, 32 x 32 grid 6712 5.3367 9607 4.8920 1.0909 0.0001 0.0001 3.7825 0.0518 True True True
overdispersion-cells (R2), mean of the per-cell variance of z, one-tailed WT > Fmr1-null:
reading n_WT WT WT_sem n_KO KO KO_sem p_t_one p_U_one WT_in_tol KO_in_tol
0 literal rule, as above 55 3.6828 0.2409 70 3.5280 0.2134 0.3157 0.2728 False False
1 literal rule, refractory filter 54 3.7093 0.2440 70 3.5280 0.2134 0.2884 0.2426 False False
2 literal rule, 32 x 32 grid 65 3.8410 0.2590 81 3.7695 0.2093 0.4142 0.4937 False False
3 count rule, as above 55 5.3285 0.4593 70 4.6888 0.3268 0.1231 0.1605 False False
4 count rule, refractory filter 54 5.3699 0.4660 70 4.6888 0.3268 0.1102 0.1421 False False
5 count rule, 32 x 32 grid 65 5.5730 0.4715 81 5.0501 0.3279 0.1755 0.2840 True True
fig4b-covariance-30-passes and -40-passes (R2), pair r over joint passes:
reading min_passes joint_min n_WT r_WT n_KO r_KO t p_t shuf1_p shuf100_p WT_in_tol KO_in_tol 0 literal rule, as above 30 39 122 0.0608 210 0.0210 2.2890 0.0227 0.2130 0.2193 True False 1 literal rule, as above 40 39 122 0.0608 209 0.0200 2.3264 0.0206 0.2279 0.2539 False False 2 literal rule, refractory filter 30 39 114 0.0606 210 0.0210 2.2777 0.0234 0.2825 0.0884 True False 3 literal rule, refractory filter 40 39 114 0.0606 209 0.0200 2.3141 0.0213 0.3003 0.1051 False False 4 literal rule, 32 x 32 grid 30 41 183 0.0744 290 0.0233 3.3426 0.0009 0.8379 0.0204 True False 5 literal rule, 32 x 32 grid 40 41 183 0.0744 290 0.0233 3.3426 0.0009 0.8379 0.0204 False False 6 count rule, as above 30 3 50 0.1646 103 0.0074 3.5426 0.0005 0.0004 0.5060 False False 7 count rule, as above 40 3 30 0.2786 86 0.0394 4.4072 0.0000 0.0639 0.9557 False False 8 count rule, refractory filter 30 3 49 0.1658 103 0.0074 3.5488 0.0005 0.0007 0.5154 False False 9 count rule, refractory filter 40 3 29 0.2846 86 0.0394 4.4669 0.0000 0.1040 0.9423 False False 10 count rule, 32 x 32 grid 30 3 94 0.1221 151 0.0264 3.2748 0.0012 0.5985 0.1294 False False 11 count rule, 32 x 32 grid 40 3 65 0.1819 121 0.0556 3.6077 0.0004 0.5340 0.4207 True True
fig4cd-kendall-place-cells, refractory filter: per reading, timescales [0.025, 0.04, 0.125] with p < 0.05 and Fmr1-null lower, and [1.0, 5.0] with p >= 0.05
K period unit short_sig long_ns n_WT n_KO
0 K1 first5 R1 0 0 110 152
1 K1 first5 R2 0 1 114 210
2 K1 session R1 0 0 110 152
3 K1 session R2 0 0 114 210
4 K2 first5 R1 0 2 106 135
5 K2 first5 R2 0 2 110 191
6 K2 session R1 0 2 110 149
7 K2 session R2 0 2 114 207
K1, session, R2: 25 ms: WT 0.0308, KO 0.0067, p = 0.845; 40 ms: WT 0.0336, KO 0.0084, p = 0.622; 125 ms: WT 0.0385, KO 0.0082, p = 0.668; 250 ms: WT 0.0416, KO 0.0079, p = 0.429; 1000 ms: WT 0.0525, KO 0.0059, p = 0.0425; 5000 ms: WT 0.0801, KO 0.0032, p = 0.00094
fig4cd-kendall-place-cells, 32 x 32 grid: per reading, timescales [0.025, 0.04, 0.125] with p < 0.05 and Fmr1-null lower, and [1.0, 5.0] with p >= 0.05
K period unit short_sig long_ns n_WT n_KO
0 K1 first5 R1 0 0 170 246
1 K1 first5 R2 0 0 183 290
2 K1 session R1 0 0 170 246
3 K1 session R2 0 0 183 290
4 K2 first5 R1 0 2 162 223
5 K2 first5 R2 0 2 175 265
6 K2 session R1 0 2 170 243
7 K2 session R2 0 1 183 287
K1, session, R2: 25 ms: WT 0.0253, KO 0.0074, p = 0.948; 40 ms: WT 0.0288, KO 0.0096, p = 0.963; 125 ms: WT 0.0346, KO 0.0108, p = 0.281; 250 ms: WT 0.0380, KO 0.0119, p = 0.242; 1000 ms: WT 0.0514, KO 0.0130, p = 0.0215; 5000 ms: WT 0.0740, KO 0.0124, p = 0.00105
kendall-nonspatial-pyramidal, refractory filter: per reading, timescales [0.025, 0.04, 0.125, 1.0] with p < 0.05 and Fmr1-null lower, and [5.0] with p >= 0.05
K period unit short_sig long_ns n_WT n_KO
0 K1 first5 R1 0 0 192 111
1 K1 first5 R2 0 0 361 132
2 K1 session R1 0 1 192 111
3 K1 session R2 1 1 361 132
4 K2 first5 R1 0 1 147 94
5 K2 first5 R2 0 1 287 114
6 K2 session R1 2 0 184 111
7 K2 session R2 2 0 350 132
K1, session, R2: 25 ms: WT 0.0092, KO 0.0101, p = 0.00902; 40 ms: WT 0.0113, KO 0.0119, p = 0.00248; 125 ms: WT 0.0157, KO 0.0152, p = 0.0112; 250 ms: WT 0.0190, KO 0.0173, p = 0.039; 1000 ms: WT 0.0242, KO 0.0196, p = 0.164; 5000 ms: WT 0.0264, KO 0.0190, p = 0.645
kendall-nonspatial-pyramidal, 32 x 32 grid: per reading, timescales [0.025, 0.04, 0.125, 1.0] with p < 0.05 and Fmr1-null lower, and [5.0] with p >= 0.05
K period unit short_sig long_ns n_WT n_KO
0 K1 first5 R1 0 1 176 72
1 K1 first5 R2 0 1 327 103
2 K1 session R1 0 1 176 72
3 K1 session R2 0 1 327 103
4 K2 first5 R1 0 1 131 57
5 K2 first5 R2 0 1 253 86
6 K2 session R1 0 0 168 72
7 K2 session R2 0 0 316 103
K1, session, R2: 25 ms: WT 0.0083, KO 0.0110, p = 0.00378; 40 ms: WT 0.0104, KO 0.0136, p = 0.000474; 125 ms: WT 0.0150, KO 0.0193, p = 0.00333; 250 ms: WT 0.0186, KO 0.0235, p = 0.00384; 1000 ms: WT 0.0242, KO 0.0301, p = 0.0178; 5000 ms: WT 0.0274, KO 0.0333, p = 0.245
Variants computed in 44 s
Under the claims file's literal pass rule neither variant changes these results. overdispersion-passes: with the refractory filter the variances are 3.71 and 3.55 (Levene p = 0.86); on the 32 x 32 grid 3.73 and 3.73 (F one-sided p = 0.49); none is within tolerance. overdispersion-cells: one-tailed p = 0.29 (filter) and 0.41 (grid), no mean within tolerance. fig4b covariance: with the filter r = 0.061 against 0.021, p = 0.023; on the grid 183 WT and 290 Fmr1-null pairs, r = 0.074 against 0.023, p = 0.0009, where the fewest joint passes is 41, so the 30 and 40 pass minimums select the same pairs; the Fmr1-null average is outside tolerance in every case, and on the grid the mean of 100 permutations still differs (p = 0.020) where one permutation per pair does not (p = 0.84). fig4cd Kendall place cells: under either variant no reading has a 25, 40 or 125 ms comparison with p < 0.05 and Fmr1-null lower (K1, session: 1 s and 5 s p = 0.043 and 0.0009 with the filter, 0.022 and 0.0011 on the grid). Kendall non-spatial pairs: at most two of the four timescales with Fmr1-null lower at p < 0.05 (K2 over the session, filter), none on the grid, where Fmr1-null tau is the higher (K1, session, 25 ms, 0.0110 against 0.0083).
Under the count rule, a variant asked of a reader: overdispersion-passes with the refractory filter has both variances within tolerance (5.30, 4.66; F one-sided p < 0.0001, Levene p = 0.0004), as without it. On the grid the variances (5.34, 4.89) and the pass counts (6,712, 9,607) are within tolerance and the F test gives p = 0.0001, but Levene's test p = 0.052. overdispersion-cells on the grid has both means within tolerance (5.57, 5.05) with one-tailed p = 0.18. fig4b-covariance-40-passes on the grid has 65 WT and 121 Fmr1-null pairs, r = 0.182 against 0.056, both within tolerance, t p = 0.0004, which meets its criterion; at 30 passes the Fmr1-null average (0.026) is outside tolerance.
kendall-quiet-states¶
During quiet wakefulness and sleep, place cell pair correlations do not differ at the gamma timescales (25, 40 ms) but are weaker in Fmr1-null at 125 ms, 250 ms and 1 s.
Criterion: tau (K1; K2 a variant) for place cell pairs over quiet periods; the readings are the recordings the data label as sleep or rest (there are none in the data set) and periods of speed below 2 cm/s in the session. Reproduced when 25 and 40 ms give p >= 0.05 and 125 ms, 250 ms and 1 s give p < 0.05 with Fmr1-null lower. Speed from room-frame positions smoothed over 1 s.
In [24]:
SPEED = {s.id: hc16.speed(s) for s in sessions}
def slow_mask(s, sample, edges):
"""Bins in which the speed is below 2 cm/s at every position sample and at the bin centre."""
v = np.nan_to_num(SPEED[s.id], nan=np.inf)
lo = np.searchsorted(s.sample_time, edges[:-1]); hi = np.searchsorted(s.sample_time, edges[1:])
vmax = np.where(hi > lo, np.maximum.reduceat(v, np.clip(lo, 0, len(v) - 1)), -np.inf)
return np.maximum(vmax, v[sample]) < 2.0
def slow_mask_centre(s, sample, edges):
return np.nan_to_num(SPEED[s.id][sample], nan=np.inf) < 2.0
out = []
for mname, mfn in [("whole bin", slow_mask), ("centre", slow_mask_centre)]:
for gate in [False, True]:
for sid, g in R2[R2.cls == "place"].groupby("session"):
cl = list(g.cell)
if len(cl) >= 2:
d = hc16.pair_taus(S[sid], cl, gate=gate, mask_fn=mfn)
d["K"] = ("K2" if gate else "K1") + ("" if mname == "whole bin" else ", centre")
out.append(d)
KQ = pd.concat(out, ignore_index=True)
KQ["in_R1"] = KQ.session.isin(set(R1.session))
slow_frac = pd.Series({sid: np.mean(np.nan_to_num(SPEED[sid], nan=np.inf) < 2) for sid in R2.session.unique()})
print("Fraction of samples below 2 cm/s per session: median %.3f, range %.3f to %.3f" % (slow_frac.median(), slow_frac.min(), slow_frac.max()))
print("Bins entering (K1, median by width), whole bin:", KQ[KQ.K == "K1"].groupby("width").nbins.median().to_dict())
print("Bins entering (K1, median by width), centre only:", KQ[KQ.K == "K1, centre"].groupby("width").nbins.median().to_dict())
rows = []
for K in ["K1", "K2", "K1, centre", "K2, centre"]:
for unit in ["R2", "R1"]:
q = KQ[KQ.K == K]
if unit == "R1":
q = q[q.in_R1]
for w in hc16.TIMESCALES:
qq = q[q.width == w].dropna(subset=["tau"])
a = qq[qq.genotype == "WT"].tau; b = qq[qq.genotype == "Fmr1-null"].tau
u = stats.mannwhitneyu(a, b)
rows.append(dict(K=K, unit=unit, width=w, n_WT=len(a), WT=a.mean(), n_KO=len(b), KO=b.mean(), U=u.statistic,
p=u.pvalue, KO_lower=b.mean() < a.mean()))
kq = pd.DataFrame(rows)
print(kq.round(4).to_string())
for key, g in kq.groupby(["K", "unit"]):
g = g.set_index("width")
print(key, "25/40 ms p >= 0.05:", int((g.loc[[0.025, 0.040], "p"] >= 0.05).sum()), "of 2;",
"125/250 ms/1 s p < 0.05 and Fmr1-null lower:", int(((g.loc[[0.125, 0.25, 1.0], "p"] < 0.05) & g.loc[[0.125, 0.25, 1.0], "KO_lower"]).sum()), "of 3")
Fraction of samples below 2 cm/s per session: median 0.212, range 0.127 to 0.419
Bins entering (K1, median by width), whole bin: {0.025: 9857.0, 0.04: 6053.0, 0.125: 1613.0, 0.25: 630.0, 1.0: 64.0, 5.0: 2.0}
Bins entering (K1, median by width), centre only: {0.025: 10179.0, 0.04: 6353.0, 0.125: 2033.0, 0.25: 1023.0, 1.0: 256.0, 5.0: 56.0}
K unit width n_WT WT n_KO KO U p KO_lower
0 K1 R2 0.025 122 0.0304 210 0.0138 12570.0 0.7764 True
1 K1 R2 0.040 122 0.0328 210 0.0185 12178.0 0.4539 True
2 K1 R2 0.125 122 0.0397 210 0.0249 12112.0 0.4081 True
3 K1 R2 0.250 122 0.0475 210 0.0303 12318.0 0.5600 True
4 K1 R2 1.000 122 0.0554 210 0.0238 13374.0 0.5039 True
5 K1 R2 5.000 50 0.2269 45 0.0925 1331.5 0.1245 True
6 K1 R1 0.025 118 0.0315 152 0.0096 9381.0 0.5169 True
7 K1 R1 0.040 118 0.0340 152 0.0133 9130.0 0.7997 True
8 K1 R1 0.125 118 0.0408 152 0.0184 8950.0 0.9781 True
9 K1 R1 0.250 118 0.0491 152 0.0218 9217.0 0.6962 True
10 K1 R1 1.000 118 0.0574 152 0.0119 9746.0 0.2218 True
11 K1 R1 5.000 50 0.2269 42 0.1106 1227.5 0.1650 True
12 K2 R2 0.025 114 0.0327 203 0.0161 11140.0 0.5825 True
13 K2 R2 0.040 115 0.0335 203 0.0234 10113.0 0.0478 True
14 K2 R2 0.125 113 0.0045 199 0.0198 9505.0 0.0232 False
15 K2 R2 0.250 103 0.0317 196 0.0200 9461.0 0.3733 True
16 K2 R2 1.000 58 0.0839 105 0.0766 3147.5 0.7237 True
17 K2 R2 5.000 0 NaN 1 0.3333 NaN NaN False
18 K2 R1 0.025 110 0.0341 146 0.0106 8178.0 0.8014 True
19 K2 R1 0.040 111 0.0351 146 0.0160 7460.0 0.2764 True
20 K2 R1 0.125 109 0.0041 142 0.0102 6900.0 0.1414 False
21 K2 R1 0.250 99 0.0343 139 0.0079 6833.0 0.9285 True
22 K2 R1 1.000 55 0.0955 61 0.1032 1724.0 0.7992 False
23 K2 R1 5.000 0 NaN 1 0.3333 NaN NaN False
24 K1, centre R2 0.025 122 0.0303 210 0.0135 12554.0 0.7619 True
25 K1, centre R2 0.040 122 0.0329 210 0.0179 12200.0 0.4698 True
26 K1, centre R2 0.125 122 0.0379 210 0.0229 12079.0 0.3863 True
27 K1, centre R2 0.250 122 0.0433 210 0.0236 12724.0 0.9192 True
28 K1, centre R2 1.000 122 0.0487 210 0.0247 13049.0 0.7773 True
29 K1, centre R2 5.000 122 0.0865 210 0.0085 15353.0 0.0026 True
30 K1, centre R1 0.025 118 0.0314 152 0.0092 9364.0 0.5343 True
31 K1, centre R1 0.040 118 0.0342 152 0.0129 9105.0 0.8302 True
32 K1, centre R1 0.125 118 0.0392 152 0.0167 8966.0 0.9981 True
33 K1, centre R1 0.250 118 0.0448 152 0.0166 9453.0 0.4465 True
34 K1, centre R1 1.000 118 0.0506 152 0.0122 9670.0 0.2704 True
35 K1, centre R1 5.000 118 0.0903 152 0.0087 10866.0 0.0029 True
36 K2, centre R2 0.025 114 0.0329 203 0.0160 11125.0 0.5694 True
37 K2, centre R2 0.040 115 0.0339 203 0.0232 10288.0 0.0789 True
38 K2, centre R2 0.125 117 0.0157 201 0.0183 9935.0 0.0211 False
39 K2, centre R2 0.250 115 0.0184 204 0.0059 10875.0 0.2800 True
40 K2, centre R2 1.000 107 -0.0154 193 0.0051 10050.0 0.7024 False
41 K2, centre R2 5.000 62 0.1394 135 0.0164 4733.5 0.1403 True
42 K2, centre R1 0.025 110 0.0344 146 0.0105 8156.0 0.8306 True
43 K2, centre R1 0.040 111 0.0354 146 0.0161 7541.0 0.3415 True
44 K2, centre R1 0.125 113 0.0162 144 0.0093 7238.0 0.1292 True
45 K2, centre R1 0.250 111 0.0184 146 -0.0096 7915.0 0.7508 True
46 K2, centre R1 1.000 103 -0.0171 136 -0.0063 6993.0 0.9842 False
47 K2, centre R1 5.000 59 0.1453 88 0.0133 2977.0 0.1326 True
('K1', 'R1') 25/40 ms p >= 0.05: 2 of 2; 125/250 ms/1 s p < 0.05 and Fmr1-null lower: 0 of 3
('K1', 'R2') 25/40 ms p >= 0.05: 2 of 2; 125/250 ms/1 s p < 0.05 and Fmr1-null lower: 0 of 3
('K1, centre', 'R1') 25/40 ms p >= 0.05: 2 of 2; 125/250 ms/1 s p < 0.05 and Fmr1-null lower: 0 of 3
('K1, centre', 'R2') 25/40 ms p >= 0.05: 2 of 2; 125/250 ms/1 s p < 0.05 and Fmr1-null lower: 0 of 3
('K2', 'R1') 25/40 ms p >= 0.05: 2 of 2; 125/250 ms/1 s p < 0.05 and Fmr1-null lower: 0 of 3
('K2', 'R2') 25/40 ms p >= 0.05: 1 of 2; 125/250 ms/1 s p < 0.05 and Fmr1-null lower: 0 of 3
('K2, centre', 'R1') 25/40 ms p >= 0.05: 2 of 2; 125/250 ms/1 s p < 0.05 and Fmr1-null lower: 0 of 3
('K2, centre', 'R2') 25/40 ms p >= 0.05: 2 of 2; 125/250 ms/1 s p < 0.05 and Fmr1-null lower: 0 of 3
/tmp/ipykernel_2508/112282393.py:34: SmallSampleWarning: One or more sample arguments is too small; all returned values will be NaN. See documentation for sample size requirements. u = stats.mannwhitneyu(a, b) /tmp/ipykernel_2508/112282393.py:34: SmallSampleWarning: One or more sample arguments is too small; all returned values will be NaN. See documentation for sample size requirements. u = stats.mannwhitneyu(a, b)
Result: not reproduced. No recording is labelled as sleep or rest, so only the low-speed reading runs. Bins enter when the speed is below 2 cm/s at every position sample in them (21% of a session's samples are below 2 cm/s, median). 25 and 40 ms do not differ (K1, R2: p = 0.78, 0.45), as claimed, but neither do 125 ms, 250 ms and 1 s (p = 0.41, 0.56, 0.50); 5 s keeps a few bins in few pairs (50 WT, 45 Fmr1-null pairs; p = 0.12). Keeping a bin when the speed at its centre is low gives the same pattern (p = 0.39, 0.92, 0.78 at 125 ms, 250 ms and 1 s). With K2 (R2), 40 ms differs (p = 0.048) and 125 ms differs with Fmr1-null higher in mean tau (p = 0.023).
locomotion-similar¶
WT and Fmr1-null mice do not differ in locomotion speed or path linearity during the recordings, and the change in speed between the two halves of a session does not differ between genotypes.
Criterion: a two-way genotype x half ANOVA with repeated measures on half, on average speed (Figure S5B inset); reproduced when the genotype effect and the interaction have p >= 0.05. Recordings and mice as units, both run. Recordings: the stationary sessions used for the Kendall claims. Linearity is not judged.
In [25]:
import pingouin as pg
rows = []
for s in sessions:
if s.env == "rotating_cylinder":
continue
v = SPEED[s.id]; h = s.nsamp // 2
rows += [dict(session=s.id, mouse=s.mouse, genotype=s.genotype, env=s.env, half="first", speed=np.nanmean(v[:h])),
dict(session=s.id, mouse=s.mouse, genotype=s.genotype, env=s.env, half="second", speed=np.nanmean(v[h:]))]
loc = pd.DataFrame(rows)
print("Speed (cm/s) by genotype and half, recordings:")
print(loc.groupby(["genotype", "half"]).speed.agg(["count", "mean", "sem"]).round(3).to_string())
print("Range %.2f to %.2f cm/s" % (loc.speed.min(), loc.speed.max()))
a1 = pg.mixed_anova(data=loc, dv="speed", within="half", between="genotype", subject="session")
print("\nRecordings as units:\n", a1.round(4).to_string())
lm = loc.groupby(["mouse", "genotype", "half"]).speed.mean().reset_index()
a2 = pg.mixed_anova(data=lm, dv="speed", within="half", between="genotype", subject="mouse")
print("\nMice as units:\n", a2.round(4).to_string())
print("\nBy mouse:")
print(lm.pivot_table(index=["genotype", "mouse"], columns="half", values="speed").round(2).to_string())
for env in ["big_box", "stationary_cylinder"]:
a3 = pg.mixed_anova(data=loc[loc.env == env], dv="speed", within="half", between="genotype", subject="session")
print(f"\nVariant, {env} recordings only:\n", a3[["Source", "F", "DF1", "DF2", [c for c in a3.columns if c.startswith("p")][0]]].round(4).to_string())
print("\nSpeed (cm/s) by environment, genotype and half, recordings; change = second - first:")
sp = loc.groupby(["env", "genotype", "half"]).speed.mean().unstack()
print(sp.assign(change=sp["second"] - sp["first"]).round(2).to_string())
for env in ["big_box", "stationary_cylinder"]:
le = loc[loc.env == env].groupby(["mouse", "genotype", "half"]).speed.mean().reset_index()
a4 = pg.mixed_anova(data=le, dv="speed", within="half", between="genotype", subject="mouse")
print(f"\nVariant, {env}, mice as units ({le.mouse.nunique()} mice):\n", a4[["Source", "F", "DF1", "DF2", [c for c in a4.columns if c.startswith("p")][0]]].round(4).to_string())
lr1 = loc[loc.session.isin(set(R1.session))]
a5 = pg.mixed_anova(data=lr1, dv="speed", within="half", between="genotype", subject="session")
print(f"\nVariant, one recording per mouse and day (R1, {lr1.session.nunique()} recordings):\n", a5[["Source", "F", "DF1", "DF2", [c for c in a5.columns if c.startswith("p")][0]]].round(4).to_string())
print("\nUnits each set of recordings gives (the paper's 27 error df need 29):",
{"stationary recordings": loc.session.nunique(), "one per mouse and day": lr1.session.nunique(),
"large box": loc[loc.env == "big_box"].session.nunique(), "stationary arena": loc[loc.env == "stationary_cylinder"].session.nunique(),
"with rotating sessions": loc.session.nunique() + len([s for s in sessions if s.env == "rotating_cylinder"]),
"mice": loc.mouse.nunique()})
Speed (cm/s) by genotype and half, recordings:
count mean sem
genotype half
Fmr1-null first 17 6.306 0.391
second 17 5.657 0.351
WT first 17 6.634 0.447
second 17 5.858 0.396
Range 2.84 to 11.55 cm/s
Recordings as units:
Source SS DF1 DF2 MS F p_unc np2 eps
0 genotype 1.1928 1 32 1.1928 0.2389 0.6283 0.0074 NaN
1 half 8.6252 1 32 8.6252 22.7869 0.0000 0.4159 1.0
2 Interaction 0.0686 1 32 0.0686 0.1813 0.6731 0.0056 1.0
Mice as units:
Source SS DF1 DF2 MS F p_unc np2 eps
0 genotype 3.1709 1 13 3.1709 0.6311 0.4412 0.0463 NaN
1 half 3.9371 1 13 3.9371 10.8646 0.0058 0.4553 1.0
2 Interaction 0.3873 1 13 0.3873 1.0689 0.3200 0.0760 1.0
By mouse:
half first second
genotype mouse
Fmr1-null M10 3.63 4.54
M18 6.58 5.69
M29 7.41 7.11
M33 5.04 4.10
M34 6.88 6.41
M37 6.39 5.16
M9 6.66 6.21
WT M1 5.45 5.63
M11 11.55 9.30
M16 7.03 5.98
M28 5.02 4.52
M3 6.14 3.64
M32 6.30 5.79
M35 6.15 5.44
M36 8.07 7.92
Variant, big_box recordings only:
Source F DF1 DF2 p_unc
0 genotype 2.0364 1 13 0.1771
1 half 11.5709 1 13 0.0047
2 Interaction 6.8945 1 13 0.0210
Variant, stationary_cylinder recordings only:
Source F DF1 DF2 p_unc
0 genotype 0.7473 1 17 0.3994
1 half 26.5423 1 17 0.0001
2 Interaction 9.1504 1 17 0.0076
Speed (cm/s) by environment, genotype and half, recordings; change = second - first:
half first second change
env genotype
big_box Fmr1-null 5.23 5.10 -0.14
WT 7.17 5.80 -1.37
stationary_cylinder Fmr1-null 7.06 6.05 -1.01
WT 6.15 5.91 -0.25
Variant, big_box, mice as units (15 mice):
Source F DF1 DF2 p_unc
0 genotype 2.0364 1 13 0.1771
1 half 11.5709 1 13 0.0047
2 Interaction 6.8945 1 13 0.0210
Variant, stationary_cylinder, mice as units (10 mice):
Source F DF1 DF2 p_unc
0 genotype 0.1661 1 8 0.6943
1 half 12.1696 1 8 0.0082
2 Interaction 7.6664 1 8 0.0243
Variant, one recording per mouse and day (R1, 25 recordings):
Source F DF1 DF2 p_unc
0 genotype 0.7002 1 23 0.4113
1 half 12.6757 1 23 0.0017
2 Interaction 1.7809 1 23 0.1951
Units each set of recordings gives (the paper's 27 error df need 29): {'stationary recordings': 34, 'one per mouse and day': 25, 'large box': 15, 'stationary arena': 19, 'with rotating sessions': 41, 'mice': 15}
Result: reproduced, resting on an assumption. With recordings as units (34 stationary sessions, large box and stationary arena pooled) neither the genotype (F(1,32) = 0.24, p = 0.63) nor the interaction (F = 0.18, p = 0.67) is significant, while speed falls from the first to the second half (F = 22.8, p < 0.0001; WT 6.63 to 5.86, Fmr1-null 6.31 to 5.66 cm/s), as the paper's epoch effect. With mice as units: genotype p = 0.44, interaction p = 0.32, half p = 0.0058; with one recording per mouse and day (25): interaction p = 0.20. Within one kind of environment the interaction is significant, in opposite directions: in the large box WT slows more (by 1.37 against 0.14 cm/s; p = 0.021, 15 mice), in the stationary arena Fmr1-null slows more (1.01 against 0.25 cm/s; p = 0.0076 with recordings, 0.024 with 10 mice as units). Pooling the environments is assumed, as the STAR Methods put all four under "Constant Conditions Open-Field" and Figure S5B speaks of the open field recordings. No set of recordings here gives the 29 units the paper's 27 error degrees of freedom imply. Linearity is not judged.
Rotating arena: place cells¶
The claims for Figures 5 and S6 use the place cells recorded in the rotating arena. Which session classifies them is not stated; two readings are run. A: place cells by the classification of the rotating session itself (room frame). B: place cells by the stationary session of the same day that shares the rotating session's cut file (stationary arena 1, or 2 for M29 whose first session has another cut file), as the paper classifies cells in standard recordings.
In [26]:
rot = [s for s in sessions if s.env == "rotating_cylinder"]
rot_place = {}
for s in rot:
a = cells[(cells.session == s.id)]
A = set(a[a.cls == "place"].cell)
same = [x for x in sessions if x.env == "stationary_cylinder" and x.mouse == s.mouse and x.date == s.date and x.cut == s.cut]
b = cells[cells.session == same[0].id] if same else None
B = set(b[b.cls == "place"].cell) & set(s.cells) if b is not None else set()
rot_place[s.id] = {"A": sorted(A), "B": sorted(B), "B_from": same[0].id if same else None}
print(f"{s.id} ({s.genotype}): {len(s.cells)} clusters; place cells A {len(A)}, B {len(B)} (from {rot_place[s.id]['B_from']}); both {len(A & B)}")
for r in ["A", "B"]:
print(r, {g: sum(len(rot_place[s.id][r]) for s in rot if s.genotype == g) for g in ["WT", "Fmr1-null"]}, "(paper: WT 20, Fmr1-null 31)")
# Ipos per 133 ms interval, per cluster and frame
T_I = time.time()
IPOS = {}
for s in rot:
for c in s.cells:
for fr in ["Rm", "Ar"]:
IPOS[(s.id, c, fr, "abs")] = hc16.ipos_series(s, c, fr)
IPOS[(s.id, c, fr, "signed_ln")] = hc16.ipos_series(s, c, fr, signed=True, log=np.log)
print(f"Ipos computed in {time.time() - T_I:.0f} s; intervals per session:", {s.id: len(IPOS[(s.id, s.cells[0], 'Rm', 'abs')]) for s in rot})
print("Fraction of intervals with Ipos defined (room, arena):",
{s.id: (round(np.mean(~np.isnan(IPOS[(s.id, s.cells[0], 'Rm', 'abs')])), 3), round(np.mean(~np.isnan(IPOS[(s.id, s.cells[0], 'Ar', 'abs')])), 3)) for s in rot})
M18_100914_rotating (Fmr1-null): 24 clusters; place cells A 3, B 4 (from M18_100914_stationary1); both 1
M29_112414_rotating (Fmr1-null): 17 clusters; place cells A 1, B 1 (from M29_112414_stationary2); both 1
M33_051815_rotating (Fmr1-null): 2 clusters; place cells A 0, B 0 (from M33_051815_stationary1); both 0
M37_051115_rotating (Fmr1-null): 23 clusters; place cells A 3, B 3 (from M37_051115_stationary1); both 3
M16_090414_rotating (WT): 50 clusters; place cells A 4, B 4 (from M16_090414_stationary1); both 1
M28_032515_rotating (WT): 10 clusters; place cells A 1, B 2 (from M28_032515_stationary1); both 1
M32_041715_rotating (WT): 15 clusters; place cells A 4, B 4 (from M32_041715_stationary1); both 2
A {'WT': 9, 'Fmr1-null': 7} (paper: WT 20, Fmr1-null 31)
B {'WT': 10, 'Fmr1-null': 8} (paper: WT 20, Fmr1-null 31)
Ipos computed in 1 s; intervals per session: {'M18_100914_rotating': 13498, 'M29_112414_rotating': 13493, 'M33_051815_rotating': 13500, 'M37_051115_rotating': 13500, 'M16_090414_rotating': 13497, 'M28_032515_rotating': 13446, 'M32_041715_rotating': 13496}
Fraction of intervals with Ipos defined (room, arena): {'M18_100914_rotating': (np.float64(0.951), np.float64(0.943)), 'M29_112414_rotating': (np.float64(0.957), np.float64(0.944)), 'M33_051815_rotating': (np.float64(0.982), np.float64(0.961)), 'M37_051115_rotating': (np.float64(0.967), np.float64(0.949)), 'M16_090414_rotating': (np.float64(0.912), np.float64(0.893)), 'M28_032515_rotating': (np.float64(0.876), np.float64(0.857)), 'M32_041715_rotating': (np.float64(0.947), np.float64(0.93))}
fig5b-frame-specific-ipos¶
During arena rotation, room-frame positional information (Ipos) is greater in Fmr1-null than in WT place cells, and arena-frame Ipos only marginally greater.
Criterion: |Ipos| per 133 ms interval from room and arena positions (32 x 32 maps, pixels with fewer than 10 visits left out), averaged over time per cell; Mann-Whitney between genotypes; reproduced when the room-frame comparison gives p < 0.05 with Fmr1-null higher. The code's signed natural-log values are a variant.
In [27]:
rows = []
for s in rot:
for c in s.cells:
rows.append(dict(session=s.id, genotype=s.genotype, mouse=s.mouse, cell=c,
A=c in rot_place[s.id]["A"], B=c in rot_place[s.id]["B"],
room=np.nanmean(IPOS[(s.id, c, "Rm", "abs")]), arena=np.nanmean(IPOS[(s.id, c, "Ar", "abs")]),
room_signed_ln=np.nanmean(IPOS[(s.id, c, "Rm", "signed_ln")]),
arena_signed_ln=np.nanmean(IPOS[(s.id, c, "Ar", "signed_ln")])))
ip = pd.DataFrame(rows)
rows = []
for r in ["A", "B"]:
d = ip[ip[r]]
for m in ["room", "arena", "room_signed_ln", "arena_signed_ln"]:
a = d[d.genotype == "WT"][m]; b = d[d.genotype == "Fmr1-null"][m]
u = stats.mannwhitneyu(a, b)
rows.append(dict(reading=r, measure=m, n_WT=len(a), WT=a.mean(), WT_sem=a.sem(), n_KO=len(b), KO=b.mean(),
KO_sem=b.sem(), U=u.statistic, p=u.pvalue, KO_higher=b.mean() > a.mean(), d=hc16.cohen_d(b, a)))
f5b = pd.DataFrame(rows)
print(f5b.round(4).to_string())
print("\nBy mouse (reading A): place cells, mean room and arena Ipos")
print(ip[ip.A].groupby(["genotype", "mouse"])[["room", "arena"]].agg(["count", "mean"]).round(4).to_string())
print("Room minus arena Ipos, place cells A: WT %.4f, Fmr1-null %.4f; cells with room > arena: %d of %d"
% ((ip[ip.A & (ip.genotype == 'WT')].room - ip[ip.A & (ip.genotype == 'WT')].arena).mean(),
(ip[ip.A & (ip.genotype == 'Fmr1-null')].room - ip[ip.A & (ip.genotype == 'Fmr1-null')].arena).mean(),
int((ip[ip.A].room > ip[ip.A].arena).sum()), int(ip.A.sum())))
print("\nPer cell, mean |Ipos| (bits) in the room and arena frame, place cells of reading A or B:")
print(ip[ip.A | ip.B][["genotype", "mouse", "session", "cell", "A", "B", "room", "arena"]].sort_values(["genotype", "room"]).round(4).to_string())
print("\nReading B, room frame, leaving out one mouse at a time: Mann-Whitney p")
d = ip[ip.B]
print({m: round(stats.mannwhitneyu(d[(d.mouse != m) & (d.genotype == "WT")].room, d[(d.mouse != m) & (d.genotype == "Fmr1-null")].room).pvalue, 4)
for m in sorted(d.mouse.unique())})
reading measure n_WT WT WT_sem n_KO KO KO_sem U p KO_higher d
0 A room 9 0.0841 0.0092 7 0.1030 0.0160 20.0 0.2523 True 0.5470
1 A arena 9 0.0809 0.0096 7 0.0920 0.0131 25.0 0.5360 True 0.3523
2 A room_signed_ln 9 0.0185 0.0031 7 0.0262 0.0062 24.0 0.4698 True 0.5963
3 A arena_signed_ln 9 0.0179 0.0035 7 0.0215 0.0044 25.0 0.5360 True 0.3209
4 B room 10 0.0663 0.0094 8 0.1056 0.0133 12.0 0.0117 True 1.1760
5 B arena 10 0.0631 0.0089 8 0.0966 0.0110 13.0 0.0155 True 1.1345
6 B room_signed_ln 10 0.0133 0.0028 8 0.0259 0.0056 17.0 0.0434 True 1.0233
7 B arena_signed_ln 10 0.0123 0.0030 8 0.0223 0.0041 20.0 0.0831 True 0.9611
By mouse (reading A): place cells, mean room and arena Ipos
room arena
count mean count mean
genotype mouse
Fmr1-null M18 3 0.0626 3 0.0590
M29 1 0.1277 1 0.1262
M37 3 0.1353 3 0.1135
WT M16 4 0.0850 4 0.0842
M28 1 0.1122 1 0.1154
M32 4 0.0761 4 0.0689
Room minus arena Ipos, place cells A: WT 0.0032, Fmr1-null 0.0110; cells with room > arena: 14 of 16
Per cell, mean |Ipos| (bits) in the room and arena frame, place cells of reading A or B:
genotype mouse session cell A B room arena
21 Fmr1-null M18 M18_100914_rotating P9C12 True False 0.0446 0.0437
17 Fmr1-null M18 M18_100914_rotating P9C8 False True 0.0550 0.0549
20 Fmr1-null M18 M18_100914_rotating P9C11 False True 0.0685 0.0646
19 Fmr1-null M18 M18_100914_rotating P9C10 True False 0.0715 0.0643
18 Fmr1-null M18 M18_100914_rotating P9C9 True True 0.0716 0.0691
54 Fmr1-null M37 M37_051115_rotating P3C2 True True 0.1147 0.0924
14 Fmr1-null M18 M18_100914_rotating P9C5 False True 0.1162 0.1171
58 Fmr1-null M37 M37_051115_rotating P3C6 True True 0.1234 0.1102
33 Fmr1-null M29 M29_112414_rotating P9C3 True True 0.1277 0.1262
44 Fmr1-null M37 M37_051115_rotating P0C2 True True 0.1678 0.1379
130 WT M32 M32_041715_rotating P1C5 False True 0.0270 0.0254
100 WT M16 M16_090414_rotating P8C9 False True 0.0404 0.0410
75 WT M16 M16_090414_rotating P0C10 True True 0.0483 0.0459
117 WT M28 M28_032515_rotating P0C2 False True 0.0501 0.0503
136 WT M32 M32_041715_rotating P1C11 False True 0.0564 0.0531
98 WT M16 M16_090414_rotating P8C7 False True 0.0580 0.0536
73 WT M16 M16_090414_rotating P0C8 True False 0.0614 0.0603
128 WT M32 M32_041715_rotating P1C3 True False 0.0615 0.0520
129 WT M32 M32_041715_rotating P1C4 True True 0.0624 0.0596
126 WT M32 M32_041715_rotating P1C1 True False 0.0756 0.0735
97 WT M16 M16_090414_rotating P8C6 False True 0.1030 0.0963
133 WT M32 M32_041715_rotating P1C8 True True 0.1050 0.0907
92 WT M16 M16_090414_rotating P8C1 True False 0.1111 0.1158
124 WT M28 M28_032515_rotating P8C1 True True 0.1122 0.1154
114 WT M16 M16_090414_rotating P11C1 True False 0.1193 0.1147
Reading B, room frame, leaving out one mouse at a time: Mann-Whitney p
{'M16': np.float64(0.0426), 'M18': np.float64(0.002), 'M28': np.float64(0.0148), 'M29': np.float64(0.025), 'M32': np.float64(0.0293), 'M37': np.float64(0.1292)}
Result: partially reproduced. Room-frame Ipos is higher in Fmr1-null under both readings, but the difference is significant in one, and rests on one mouse. With place cells classified in the stationary session of the same day (reading B; 10 WT, 8 Fmr1-null) room-frame Ipos is greater in Fmr1-null (0.106 against 0.066 bits, U = 12, p = 0.012); arena-frame Ipos is greater too (0.097 against 0.063, p = 0.016), where the paper finds it only marginal. Leaving out M37 (3 of the 8 Fmr1-null cells) gives p = 0.13; leaving out any other mouse keeps p below 0.05. With place cells classified in the rotating session itself (reading A; 9 and 7) the room-frame difference has the same direction but is not significant (0.103 against 0.084, p = 0.25). The samples come from 3 mice per genotype and are about half the paper's 20 and 31 cells.
fig5c-ensemble-coordination¶
The frame-specific positional discharge fluctuations of single place cells are less correlated with those of the rest of the ensemble in Fmr1-null than in WT mice.
Criterion: per cell, delta Ipos = Ipos(room) - Ipos(arena) per 133 ms interval; the ensemble series from the frame-specific sums of Ipos over the other simultaneously recorded cells; the Pearson correlation of the cell's and the ensemble's series; delta and |delta| both run; t test on Fisher z; reproduced when, under one reading, p < 0.05 with Fmr1-null lower. The ensemble is "all simultaneously recorded cells except one" (Results); the other place cells only is a variant. Place cell readings A and B.
In [28]:
rows = []
for s in rot:
allc = s.cells
room = np.vstack([IPOS[(s.id, c, "Rm", "abs")] for c in allc])
arena = np.vstack([IPOS[(s.id, c, "Ar", "abs")] for c in allc])
for r in ["A", "B"]:
pcs = rot_place[s.id][r]
for c in pcs:
i = allc.index(c)
for ens_name, members in [("all clusters", [k for k in range(len(allc)) if k != i]),
("place cells", [allc.index(x) for x in pcs if x != c])]:
if len(members) == 0:
continue
d_cell = room[i] - arena[i]
d_ens = np.nansum(room[members], 0) - np.nansum(arena[members], 0)
ok = ~np.isnan(d_cell) & ~np.isnan(room[members]).all(0) & ~np.isnan(arena[members]).all(0)
for form in ["delta", "abs delta"]:
x = d_cell[ok]; y = d_ens[ok]
if form == "abs delta":
x, y = np.abs(x), np.abs(y)
rr = stats.pearsonr(x, y).statistic if x.std() > 0 and y.std() > 0 else np.nan
rows.append(dict(session=s.id, genotype=s.genotype, mouse=s.mouse, cell=c, reading=r, ensemble=ens_name,
form=form, n_ens=len(members), r=rr))
ec = pd.DataFrame(rows)
rows = []
for key, g in ec.groupby(["reading", "ensemble", "form"]):
a = g[g.genotype == "WT"].r.dropna(); b = g[g.genotype == "Fmr1-null"].r.dropna()
if len(a) < 2 or len(b) < 2:
continue
t = stats.ttest_ind(hc16.fisher(a), hc16.fisher(b))
rows.append(dict(reading=key[0], ensemble=key[1], form=key[2], n_WT=len(a), r_WT=a.mean(), sem_WT=a.sem(),
n_KO=len(b), r_KO=b.mean(), sem_KO=b.sem(), t=t.statistic, p=t.pvalue, KO_lower=b.mean() < a.mean(),
d=hc16.cohen_d(hc16.fisher(a), hc16.fisher(b))))
f5c = pd.DataFrame(rows)
print(f5c.round(4).to_string())
print("\nBy mouse (A, all clusters, abs delta):")
print(ec[(ec.reading == "A") & (ec.ensemble == "all clusters") & (ec.form == "abs delta")].groupby(["genotype", "mouse"]).r.agg(["count", "mean"]).round(3).to_string())
print("r range: %.3f to %.3f" % (ec.r.min(), ec.r.max()))
reading ensemble form n_WT r_WT sem_WT n_KO r_KO sem_KO t p KO_lower d
0 A all clusters abs delta 9 0.0704 0.0214 7 0.0713 0.0329 -0.0289 0.9773 False -0.0146
1 A all clusters delta 9 0.1301 0.0268 7 0.1566 0.0185 -0.7524 0.4643 False -0.3792
2 A place cells abs delta 8 0.0093 0.0132 6 0.0278 0.0301 -0.6208 0.5463 False -0.3353
3 A place cells delta 8 0.0532 0.0162 6 0.0985 0.0232 -1.6521 0.1244 False -0.8922
4 B all clusters abs delta 10 0.0447 0.0151 8 0.0759 0.0246 -1.1333 0.2738 False -0.5376
5 B all clusters delta 10 0.1147 0.0195 8 0.1752 0.0188 -2.1802 0.0445 False -1.0342
6 B place cells abs delta 10 0.0139 0.0084 7 0.0285 0.0183 -0.8052 0.4333 False -0.3968
7 B place cells delta 10 0.0451 0.0088 7 0.0919 0.0109 -3.3733 0.0042 False -1.6624
By mouse (A, all clusters, abs delta):
count mean
genotype mouse
Fmr1-null M18 3 0.086
M29 1 0.012
M37 3 0.076
WT M16 4 0.103
M28 1 0.046
M32 4 0.044
r range: -0.049 to 0.260
Result: not reproduced. Under every reading the Fmr1-null place cells correlate with the ensemble at least as well as the WT cells. With |delta Ipos| and all other clusters as the ensemble, r = 0.071 (Fmr1-null) against 0.070 (WT) in reading A (p = 0.98) and 0.076 against 0.045 in reading B (p = 0.27); with delta Ipos, Fmr1-null is higher in reading B with all clusters (0.175 against 0.115, p = 0.045) and with the other place cells (0.092 against 0.045, p = 0.0042).
figs6-room-frame-preference¶
Session-averaged spatial firing during rotation is better organized in the room frame than in the arena frame, similarly in both genotypes.
Criterion: a two-way genotype x environment ANOVA with repeated measures on environment (Stationary I, Stationary II, rotating room frame, rotating arena frame) on coherence and information content (log transformed; the reverse as a variant), with Tukey post hoc tests; reproduced when the environment effect has p < 0.05 for coherence with the rotating arena frame below both stationary sessions, and the interaction has p >= 0.05. Cells must be present in all four maps, so they must share the cut file across the day's three sessions (M29, whose stationary arena 1 has another cut file, drops out). Tukey's HSD is applied across the four environments without the repeated structure (statsmodels); paired t tests are shown as a variant.
In [29]:
import pingouin as pg
from statsmodels.stats.multicomp import pairwise_tukeyhsd
def map_measures(s, c, frame):
idx, shape = hc16.pixel_index(s, frame)
rm, o = hc16.rate_map(s, c, idx, shape)
return hc16.coherence(rm), hc16.information(rm, o)
rows = []
for s in rot:
st = {x.cond: x for x in sessions if x.env == "stationary_cylinder" and x.mouse == s.mouse and x.date == s.date and x.cut == s.cut}
if not {"stationary1", "stationary2"} <= set(st):
print(f"{s.id}: no stationary arena 1 and 2 with the same cut file; left out")
continue
for r in ["A", "B"]:
for c in rot_place[s.id][r]:
if c not in st["stationary1"].cells or c not in st["stationary2"].cells:
continue
for env, (ss, fr) in {"Stationary I": (st["stationary1"], "Rm"), "Stationary II": (st["stationary2"], "Rm"),
"Rotating room": (s, "Rm"), "Rotating arena": (s, "Ar")}.items():
coh, inf = map_measures(ss, c, fr)
rows.append(dict(reading=r, session=s.id, genotype=s.genotype, mouse=s.mouse, cell=f"{s.id}.{c}", env=env,
coherence=coh, info=inf, log_info=np.log(inf), log_coherence=np.log(coh) if coh > 0 else np.nan))
fs6 = pd.DataFrame(rows)
for r in ["A", "B"]:
d = fs6[fs6.reading == r].dropna(subset=["coherence", "log_info"])
full = d.groupby("cell").env.nunique() == 4
d = d[d.cell.isin(full[full].index)]
print(f"\nReading {r}: cells with all four maps: {d.groupby('genotype').cell.nunique().to_dict()} (paper: WT 20, Fmr1-null 31)")
print(d.groupby(["env", "genotype"])[["coherence", "info"]].mean().round(3).unstack().to_string())
for m in ["coherence", "log_info", "info", "log_coherence"]:
dd = d.dropna(subset=[m])
keep = dd.groupby("cell").env.nunique() == 4
dd = dd[dd.cell.isin(keep[keep].index)]
a = pg.mixed_anova(data=dd, dv=m, within="env", between="genotype", subject="cell")
pcol = [c for c in a.columns if c.startswith("p")][0]
print(f" {m}: " + "; ".join(f"{row.Source} F({row.DF1},{row.DF2}) = {row.F:.2f}, p = {row[pcol]:.4g}" for _, row in a.iterrows()))
if m in ["coherence", "log_info"]:
tk = pairwise_tukeyhsd(dd[m], dd.env)
print(" Tukey HSD:\n" + "\n".join(" " + l for l in str(tk.summary()).splitlines()[2:]))
w = dd.pivot(index="cell", columns="env", values=m)
for e in ["Stationary I", "Stationary II"]:
tt = stats.ttest_rel(w["Rotating arena"], w[e])
print(f" paired t, Rotating arena vs {e}: mean diff {np.mean(w['Rotating arena'] - w[e]):.3f}, t = {tt.statistic:.2f}, p = {tt.pvalue:.4g}")
M29_112414_rotating: no stationary arena 1 and 2 with the same cut file; left out
Reading A: cells with all four maps: {'Fmr1-null': 6, 'WT': 9} (paper: WT 20, Fmr1-null 31)
coherence info
genotype Fmr1-null WT Fmr1-null WT
env
Rotating arena 0.394 0.423 1.373 1.371
Rotating room 0.753 0.532 1.652 1.481
Stationary I 0.536 0.464 2.590 1.958
Stationary II 0.687 0.466 1.592 1.730
coherence: genotype F(1,13) = 2.28, p = 0.1551; env F(3,39) = 2.58, p = 0.06716; Interaction F(3,39) = 1.18, p = 0.3309
Tukey HSD:
group1 group2 meandiff p-adj lower upper reject
------------------------------------------------------------------
Rotating arena Rotating room 0.209 0.1031 -0.0284 0.4463 False
Rotating arena Stationary I 0.0816 0.7993 -0.1557 0.3189 False
Rotating arena Stationary II 0.1428 0.3909 -0.0946 0.3801 False
Rotating room Stationary I -0.1274 0.4917 -0.3647 0.11 False
Rotating room Stationary II -0.0662 0.881 -0.3035 0.1711 False
Stationary I Stationary II 0.0612 0.9034 -0.1762 0.2985 False
------------------------------------------------------------------
paired t, Rotating arena vs Stationary I: mean diff -0.082, t = -0.96, p = 0.352
paired t, Rotating arena vs Stationary II: mean diff -0.143, t = -2.10, p = 0.05396
log_info: genotype F(1,13) = 0.02, p = 0.8791; env F(3,39) = 12.68, p = 6.261e-06; Interaction F(3,39) = 2.38, p = 0.0839
Tukey HSD:
group1 group2 meandiff p-adj lower upper reject
------------------------------------------------------------------
Rotating arena Rotating room 0.1238 0.9199 -0.3913 0.6388 False
Rotating arena Stationary I 0.4359 0.1248 -0.0791 0.9509 False
Rotating arena Stationary II 0.2277 0.6478 -0.2874 0.7427 False
Rotating room Stationary I 0.3122 0.3843 -0.2029 0.8272 False
Rotating room Stationary II 0.1039 0.9503 -0.4111 0.6189 False
Stationary I Stationary II -0.2082 0.7087 -0.7233 0.3068 False
------------------------------------------------------------------
paired t, Rotating arena vs Stationary I: mean diff -0.436, t = -5.16, p = 0.0001439
paired t, Rotating arena vs Stationary II: mean diff -0.228, t = -3.14, p = 0.007253
info: genotype F(1,13) = 0.16, p = 0.699; env F(3,39) = 7.03, p = 0.0006865; Interaction F(3,39) = 1.45, p = 0.2436
log_coherence: genotype F(1,13) = 0.14, p = 0.7111; env F(3,39) = 2.48, p = 0.07547; Interaction F(3,39) = 0.64, p = 0.5966
Reading B: cells with all four maps: {'Fmr1-null': 7, 'WT': 10} (paper: WT 20, Fmr1-null 31)
coherence info
genotype Fmr1-null WT Fmr1-null WT
env
Rotating arena 0.334 0.260 1.022 1.741
Rotating room 0.527 0.338 1.182 1.883
Stationary I 0.592 0.581 1.544 1.948
Stationary II 0.506 0.352 1.196 2.157
coherence: genotype F(1,15) = 1.14, p = 0.3028; env F(3,45) = 6.24, p = 0.001241; Interaction F(3,45) = 0.66, p = 0.5798
Tukey HSD:
group1 group2 meandiff p-adj lower upper reject
------------------------------------------------------------------
Rotating arena Rotating room 0.1251 0.5249 -0.1164 0.3666 False
Rotating arena Stationary I 0.2947 0.0107 0.0531 0.5362 True
Rotating arena Stationary II 0.1243 0.5301 -0.1172 0.3659 False
Rotating room Stationary I 0.1696 0.2591 -0.072 0.4111 False
Rotating room Stationary II -0.0008 1.0 -0.2423 0.2408 False
Stationary I Stationary II -0.1703 0.2555 -0.4119 0.0712 False
------------------------------------------------------------------
paired t, Rotating arena vs Stationary I: mean diff -0.295, t = -3.81, p = 0.001542
paired t, Rotating arena vs Stationary II: mean diff -0.124, t = -2.06, p = 0.05643
log_info: genotype F(1,15) = 3.00, p = 0.104; env F(3,45) = 6.60, p = 0.0008615; Interaction F(3,45) = 1.99, p = 0.1285
Tukey HSD:
group1 group2 meandiff p-adj lower upper reject
------------------------------------------------------------------
Rotating arena Rotating room 0.1087 0.95 -0.4268 0.6443 False
Rotating arena Stationary I 0.2551 0.5936 -0.2804 0.7907 False
Rotating arena Stationary II 0.2066 0.7397 -0.3289 0.7422 False
Rotating room Stationary I 0.1464 0.8884 -0.3892 0.682 False
Rotating room Stationary II 0.0979 0.9628 -0.4377 0.6335 False
Stationary I Stationary II -0.0485 0.9952 -0.5841 0.4871 False
------------------------------------------------------------------
paired t, Rotating arena vs Stationary I: mean diff -0.255, t = -3.27, p = 0.004766
paired t, Rotating arena vs Stationary II: mean diff -0.207, t = -3.21, p = 0.00549
info: genotype F(1,15) = 2.79, p = 0.1154; env F(3,45) = 3.78, p = 0.01684; Interaction F(3,45) = 1.91, p = 0.1414
log_coherence: genotype F(1,13) = 0.87, p = 0.3673; env F(3,39) = 5.47, p = 0.003089; Interaction F(3,39) = 0.44, p = 0.7257
Result: partially reproduced, resting on an assumption. In reading B (10 WT and 7 Fmr1-null place cells with all four maps) coherence has an environment effect (F(3,45) = 6.24, p = 0.0012) and no interaction (p = 0.58); Tukey's test puts the rotating arena frame below stationary arena I (difference 0.29, adjusted p = 0.011) but not below stationary arena II (0.12, p = 0.53; paired t p = 0.056). Log information has an environment effect (p = 0.0009) and no interaction (p = 0.13). In reading A (9 and 6 cells) the coherence environment effect is not significant (p = 0.067). The outcome adopts reading B, an assumption. It settles a question the record leaves open elsewhere: fig5b-frame-specific-ipos is judged on both readings, and the question is put to the paper's authors. Under reading A this claim would not be reproduced.
Summary¶
| Claim | Outcome |
|---|---|
| fig1ab-baseline-transmission | blocked: data missing (slice recordings) |
| fig1cd-potentiation | blocked: data missing (slice recordings) |
| figs1-avoidance-equal | blocked: data missing (slice-cohort behaviour) |
| unit-and-cell-counts | blocked: data missing (15 of 21 mice; no isolation measures) |
| interneuron-proportion | partially reproduced (direction holds; percentages far higher) |
| place-cell-proportion | not reproduced |
| figs3-nonspatial-properties-similar | partially reproduced (rate differs) |
| fig2-spatial-quality-similar | not reproduced (genotype effects, mostly in interneurons) |
| tables1-place-cell-properties-similar | partially reproduced (place cells alike; interneuron rate differs) |
| tables1-nonspatial-active-pixels | not reproduced |
| figs4-remapping-similar | blocked: data missing (no square and circular recordings of the same cells) |
| fig3b-place-cell-modulation-similar, fig3bc-weaker-modulation-nonspatial-interneurons, fig3b-phase-preference, fig3d-stereotypy, figs5a-power-spectra-similar | blocked: data missing (no LFP) |
| overdispersion-passes | not reproduced (literal pass rule; the count rule meets the criterion in R2, asked of a reader) |
| overdispersion-cells | not reproduced |
| fig4b-covariance-30-passes | partially reproduced |
| fig4b-covariance-40-passes | partially reproduced |
| fig4cd-kendall-place-cells | not reproduced |
| fig4d-network-state-recurrence | not reproduced |
| kendall-nonspatial-pyramidal | not reproduced |
| kendall-quiet-states | not reproduced |
| locomotion-similar | reproduced (assumption: environments pooled) |
| fig5b-frame-specific-ipos | partially reproduced (significant in one reading, resting on one mouse) |
| fig5c-ensemble-coordination | not reproduced |
| figs6-room-frame-preference | partially reproduced (assumption: cells classified in the stationary session) |
| fig6-covariance-sem-predicts-errors, fig6-kendall-predicts-errors, discussion-single-cell-measures-no-relation, fig6-avoidance-behavior | blocked: data missing (avoidance sessions) |
A failure to reproduce here may be a failure of this attempt, not of the paper: see the reasoning and the differences from the paper's method.
How the analysis differs from the paper's:
- The paper's cells are units with IsoIBG and IsoINN of at least 4 bits. The data set has no isolation measures and no waveforms, so every sorted cluster other than cluster 0 is used, and cells are classified by rate and bursting without waveform width. Table S2 lists fewer cells per session than the data set has clusters (for example 15 cells for M16 in stationary arena 1, against 50 clusters here), so the population here includes clusters the paper excluded. This bears on every tested claim. A refractory filter, a weak proxy for isolation, removes 70 of 774 clusters. Among the single-cell claims it removes the rate effect of figs3-nonspatial-properties-similar and, in R1, the significance of the interneuron difference, and leaves the other results as they are. Run on the pass and Kendall claims (R2), it leaves their results as they are. It was not run on fig4d-network-state-recurrence, kendall-quiet-states, locomotion-similar or the rotating arena claims.
- The data set holds 8 WT and 7 Fmr1-null mice, the large box and the stationary and rotating arena. The small box, the ring, sleep recordings, LFP, slice data and avoidance sessions are not in it.
- The spatial scale is assumed from the tracked extent and the enclosure sizes of the Methods. On the code's 32 x 32 grid, which needs no scale, 45 of 774 clusters change class, and the outcome of no single-cell claim changes. On the grid, under the literal pass rule, the results of the pass and Kendall claims (R2) are as on 1.25 cm pixels; under the count rule fig4b-covariance-40-passes meets its criterion. The grid was not run on fig4d-network-state-recurrence, kendall-quiet-states, locomotion-similar or the rotating arena claims (there Ipos is computed on the 32 x 32 grid, but the place cells are classified on 1.25 cm pixels).
- The pass rule of the overdispersion claims is read literally, as the claims file does; it passes most 5 s intervals. The count expected at the mean rate over 5 s, a variant, changes the overdispersion-passes result in R2 (not in R1, and not overdispersion-cells), and with the 32 x 32 grid the fig4b-covariance-40-passes result; it is asked of a reader.
- The session that classifies a cell is open (R1, R2; A, B for the rotating arena); the results under each are printed. A cluster name within a day's cut file is taken to be one cell. R2 counts the cells and pairs of stationary arena 1 and 2 of a day twice; R1 counts them once.
What was found about the data set:
- It holds part of the paper's data: 15 of 21 mice, three of the four environments, no LFP, no sleep or rest recordings, no avoidance sessions, and no slice recordings.
- No isolation measures, waveforms or cell classes are given, so the paper's cells cannot be selected; the files hold every cluster, including cluster 0 (unsorted spikes). The TS headers hold only the sampling rate, the first timestamp and a track-to-map scale of 1.00.
- The descriptor gives no spatial scale for the 0-255 pixel positions, and the scale differs between box sessions (about 3.6 and 4.4 px/cm by extent).
- The stated position sampling rate drifts against the spike clock (spikes up to 118 samples, about 4 s, after the nominal time of their position sample); spikes are filed under the right sample, so the files are consistent, but the nominal sample times are not usable as times.
- Six spike trains have one backward step in their times (up to 1.4 s).
- The 0-15 and 15-30 min files repeat the full sessions except one position sample at the split,
which is in neither part, so a spike filed under it is lost; the
.timesfiles repeat the TS spike times. - The descriptor's example header gives 58.3 samples/s; every file has 30.0.
In [30]:
print(f"Notebook run time: {(time.time() - T0) / 60:.1f} min")
Notebook run time: 2.3 min