import time
import warnings

import matplotlib.pyplot as plt
import numpy as np
import pandas as pd

import clark2024 as c

pd.set_option("display.width", 220)
pd.set_option("display.max_columns", 40)
pd.set_option("display.max_rows", 400)
warnings.filterwarnings("ignore")

started = time.time()
info = c.version_info()
assets = c.list_assets()
print(f"dandiset {c.DANDISET_ID}, version {c.DANDISET_VERSION}: {len(assets)} assets, "
      f"{sum(a['size'] for a in assets) / 1e9:.1f} GB")
for v in info["versions"]:
    print(f"  version {v['version']:>14}: {v['asset_count']} assets, {v['size'] / 1e9:.1f} GB, "
          f"status {v['status']}, modified {v['modified'][:19]}")
meta = info["meta"]
print("\nDescription:", meta["description"])
print("variableMeasured:", [x["value"] if isinstance(x, dict) else x for x in meta["assetsSummary"].get("variableMeasured", [])])
print("numberOfSubjects:", meta["assetsSummary"].get("numberOfSubjects"))

dandiset 000943, version 0.240614.0754: 734 assets, 976.2 GB
  version          draft: 734 assets, 976.2 GB, status Published, modified 2024-06-14T07:54:16
  version  0.240614.0754: 734 assets, 976.2 GB, status Valid, modified 2024-06-14T07:54:17

Description: This dataset contains electrophysiological recordings alongside the synchronised behavioural data for mice running in an open field and virtual linear track environments. Movable tetrodes were targeted to the medial entorhinal cortex. Further details can be found here (https://doi.org/10.7554/eLife.89356.2).
variableMeasured: ['ElectrodeGroup', 'ElectricalSeries']
numberOfSubjects: 12

# One streamed open per file: metadata, the contents of every group that could hold derived data,
# and a summary of the position series (read whole; it is small).
t = time.time()
survey = pd.DataFrame(c.map_files(c.survey_file, assets, workers=16))
print(f"surveyed {len(survey)} files in {time.time() - t:.0f} s; files that failed to open: "
      f"{int(survey['error'].notna().sum()) if 'error' in survey else 0}")
survey["day"] = survey.session_id.str.split("_").str[1]
vr = survey[survey.session_description == "vr"].copy()
of = survey[survey.session_description == "openfield"].copy()

surveyed 734 files in 83 s; files that failed to open: 0

print("Files by mouse and session type (session_description):")
print(pd.crosstab(survey.mouse, survey.session_description, margins=True).to_string())

print("\nWhat the files hold (every file):")
contents = survey.contents.map(lambda d: str({k: v for k, v in d.items()}))
print(contents.value_counts().to_string())
n_units = sum(d["units"] is not None for d in survey.contents)
n_proc = sum(bool(d["processing"]) for d in survey.contents)
n_intervals = sum(d["intervals"] is not None for d in survey.contents)
print(f"files with a units table: {n_units}; with any processing module: {n_proc}; with intervals (trials): {n_intervals}")

print("\nRaw signal: channels", sorted(survey.es_channels.unique()), "; rate (Hz)", sorted(survey.es_rate.unique()),
      "; conversion", sorted(survey.es_conversion.unique()), sorted(survey.es_unit.unique()))
dur = survey.es_samples / survey.es_rate / 60
print("Session length from the raw signal (min) by type:")
print(dur.groupby(survey.session_description).describe().round(1).to_string())
print("Position series description (identical in every file):", survey.pos_description.unique().tolist())
print("Position columns by type:", survey.groupby("session_description").pos_shape.agg(lambda s: sorted({p[1] for p in s})).to_dict())

Files by mouse and session type (session_description):
session_description  openfield   vr  All
mouse                                   
M1                          41   11   52
M10                         25   25   50
M11                         45   45   90
M12                         30   30   60
M13                         32   29   61
M14                         45   45   90
M15                         33   21   54
M2                          21   29   50
M3                          29   35   64
M4                          29   30   59
M6                          22   22   44
M7                          27   33   60
All                        379  355  734

What the files hold (every file):
contents
{'processing': [], 'analysis': [], 'stimulus': ['presentation', 'templates'], 'intervals': None, 'units': None, 'scratch': None, 'acquisition': ['ElectricalSeries', 'position']}    734
files with a units table: 0; with any processing module: 0; with intervals (trials): 0

Raw signal: channels [np.int64(16)] ; rate (Hz) [np.float64(30000.0)] ; conversion [np.float64(1.9499999999999999e-07)] ['volts']
Session length from the raw signal (min) by type:
                     count  mean  std   min   25%   50%   75%   max
session_description                                                
openfield            379.0  30.4  2.1  23.3  29.7  30.2  30.9  52.9
vr                   355.0  30.6  3.2  11.3  29.7  30.2  31.1  42.6
Position series description (identical in every file): ['position: cm, speed: cm/s, angle: degrees, trials: non-zero indexed, trial types: 0-beaconed, 1-nonbeaconed, 2-probe']
Position columns by type: {'openfield': [9], 'vr': [3]}

# The free-text metadata of one VR file and one open arena file, which the checks below refer to.
import textwrap
by_path = {a["path"]: a for a in assets}
example = [vr.sort_values("path").path.iloc[0], of.sort_values("path").path.iloc[0]]
texts = c.map_files(c.file_text, [by_path[p] for p in example], workers=2)
for t in texts:
    print(f"== {t['path']}")
    for k, v in t.items():
        if k != "path":
            print(textwrap.fill(f"{k}: {v}", 120, subsequent_indent="    "))
print("\ngeneral/notes identical in the two files:", texts[0]["general/notes"] == texts[1]["general/notes"],
      "; experiment_description identical:", texts[0]["general/experiment_description"] == texts[1]["general/experiment_description"])

== sub-M1/sub-M1_ses-M1-D1-2020-08-03-16-11-14_ecephys.nwb
session_description: vr
general/notes: Experimental days involved recording from mice in the open arena and then in the virtual location memory
    task. On a few days this order was reversed without apparent effects on the results obtained. Mice were collected
    from the holding room 30 - 60 minutes before recording, were handled for 5 - 10 minutes, weighed and placed for 10 -
    20 minutes in a cage containing objects and a running wheel. Between recording sessions mice were placed back in the
    object-filled playground for 30 minutes. Tetrodes were typically lowered by 50 - 100 µm after each session. The open
    arena consisted of a metal box with a square floor area, removable metal walls, metal frame (Frame parts from Kanya
    UK, C01-1, C20-10, A33-12, B49-75, B48-75, A39-31, ALU3), and an A4-sized cue card in the middle of one of the metal
    walls. For the open field exploration session, mice were placed in the open arena while tethered via an ultrathin
    SPI cable and custom-build commutator and left unprompted for 30 minutes to freely move around. For the location
    memory task mice were trained to obtain rewards at a location on the virtual linear track. Mice were head-fixed
    using a RIVETS clamp (Ronal Tool Company, Inc) and ran on a cylindrical treadmill fitted with a rotary encoder
    (Pewatron). Virtual tracks, generated using Blender3D (blender.com) had a range of lengths, in the published journal
    article only sessions with 200cm track lengths are used. This consisted of a 60 cm track zone, a 20 cm reward zone,
    a second 60 cm track zone and a 60 cm black box to separate successive trials. In sessions with track lengths longer
    than 200cm, Xcm length tracks consisted of a Zcm length track zone, a 20 cm reward zone, a 60 cm track zone and a 60
    cm black box to seperate successive trials. The length of Z can be derived from the track length using the equations
    Z = X-(60+60+20). In actuality, the trial starts in the middle of the blackbox region therefore it can be thought of
    as 2 successive 30 cm black box regions or a single 60 cm region that spans across two trials, this explains how the
    black box is visualised within the research article. The distance visible ahead of the mouse was 50 cm. The reward
    zone was either marked by distinct vertical green and black bars on beaconed trials, or was not marked by a visual
    cue at all on non-beaconed or probe trials. A feeding tube placed in front of the animal dispensed soy milk rewards
    (5 - 10 µl per reward) if the mouse stopped in the reward zone, however was not dispensed on probe trials. A stop
    was registered in Blender3D if the speed of the mouse dropped below 4.7 cm/s. Speed was calculated on a rolling
    basis from the previous 100 ms at a rate of 60 Hz. Trials were delivered in repeating blocks throughout a recording
    session. For example, 3 beaconed trials (B) followed by 2 non-beaconed trials (N) with the order repeating until the
    end of the session. To encourage learning and engagement in both beaconed and non-beaconed trials, the first day of
    training typically used a trial type ratio of 3 beaconed trials to 1 non-beaconed trials. As training progressed we
    then increased the proportion of non-beaconed trials up to a ratio of 1 beaconed trial to 4 non-beaconed trials.
    Examples of trial blocks include BBBBN, BBBN, BBN, BN, BNN, BNNN, BBNNN and BBNNNNNNN where each character indicates
    the trial type of each trial within a block. In some sessions we replaced single non-beaconed trials in trial blocks
    with a probe trial (P). Examples of trial blocks with probe trials include BBBBNBBBBP, BBBNBBBP, BBNBBP and
    BBNNNNNNNP. For more details, see https://doi.org/10.7554/eLife.89356.2
general/experiment_description: Extracellular recording of mice medial entorhinal cortex during free exploration in a
    100 x 100 cm open field and a virtual reality location memory task
general/protocol: see https://doi.org/10.7554/eLife.89356.2
general/session_id: M1_D1_2020-08-03_16-11-14
== sub-M1/sub-M1_ses-M1-D1-2020-08-03-16-49-15_ecephys.nwb
session_description: openfield
general/notes: Experimental days involved recording from mice in the open arena and then in the virtual location memory
    task. On a few days this order was reversed without apparent effects on the results obtained. Mice were collected
    from the holding room 30 - 60 minutes before recording, were handled for 5 - 10 minutes, weighed and placed for 10 -
    20 minutes in a cage containing objects and a running wheel. Between recording sessions mice were placed back in the
    object-filled playground for 30 minutes. Tetrodes were typically lowered by 50 - 100 µm after each session. The open
    arena consisted of a metal box with a square floor area, removable metal walls, metal frame (Frame parts from Kanya
    UK, C01-1, C20-10, A33-12, B49-75, B48-75, A39-31, ALU3), and an A4-sized cue card in the middle of one of the metal
    walls. For the open field exploration session, mice were placed in the open arena while tethered via an ultrathin
    SPI cable and custom-build commutator and left unprompted for 30 minutes to freely move around. For the location
    memory task mice were trained to obtain rewards at a location on the virtual linear track. Mice were head-fixed
    using a RIVETS clamp (Ronal Tool Company, Inc) and ran on a cylindrical treadmill fitted with a rotary encoder
    (Pewatron). Virtual tracks, generated using Blender3D (blender.com) had a range of lengths, in the published journal
    article only sessions with 200cm track lengths are used. This consisted of a 60 cm track zone, a 20 cm reward zone,
    a second 60 cm track zone and a 60 cm black box to separate successive trials. In sessions with track lengths longer
    than 200cm, Xcm length tracks consisted of a Zcm length track zone, a 20 cm reward zone, a 60 cm track zone and a 60
    cm black box to seperate successive trials. The length of Z can be derived from the track length using the equations
    Z = X-(60+60+20). In actuality, the trial starts in the middle of the blackbox region therefore it can be thought of
    as 2 successive 30 cm black box regions or a single 60 cm region that spans across two trials, this explains how the
    black box is visualised within the research article. The distance visible ahead of the mouse was 50 cm. The reward
    zone was either marked by distinct vertical green and black bars on beaconed trials, or was not marked by a visual
    cue at all on non-beaconed or probe trials. A feeding tube placed in front of the animal dispensed soy milk rewards
    (5 - 10 µl per reward) if the mouse stopped in the reward zone, however was not dispensed on probe trials. A stop
    was registered in Blender3D if the speed of the mouse dropped below 4.7 cm/s. Speed was calculated on a rolling
    basis from the previous 100 ms at a rate of 60 Hz. Trials were delivered in repeating blocks throughout a recording
    session. For example, 3 beaconed trials (B) followed by 2 non-beaconed trials (N) with the order repeating until the
    end of the session. To encourage learning and engagement in both beaconed and non-beaconed trials, the first day of
    training typically used a trial type ratio of 3 beaconed trials to 1 non-beaconed trials. As training progressed we
    then increased the proportion of non-beaconed trials up to a ratio of 1 beaconed trial to 4 non-beaconed trials.
    Examples of trial blocks include BBBBN, BBBN, BBN, BN, BNN, BNNN, BBNNN and BBNNNNNNN where each character indicates
    the trial type of each trial within a block. In some sessions we replaced single non-beaconed trials in trial blocks
    with a probe trial (P). Examples of trial blocks with probe trials include BBBBNBBBBP, BBBNBBBP, BBNBBP and
    BBNNNNNNNP. For more details, see https://doi.org/10.7554/eLife.89356.2
general/experiment_description: Extracellular recording of mice medial entorhinal cortex during free exploration in a
    100 x 100 cm open field and a virtual reality location memory task
general/protocol: see https://doi.org/10.7554/eLife.89356.2
general/session_id: M1_D1_2020-08-03_16-49-15

general/notes identical in the two files: True ; experiment_description identical: True

# Virtual track sessions: columns are position (cm), trial number and trial type (description above).
print("VR position: sample interval (s)", sorted(vr.dt.round(6).unique()), "; first timestamp", sorted(vr.t_first.unique()))
print("trial numbers start at 1 in every file:", bool((vr.trial_first == 1).all()),
      "; never decrease:", bool(vr.trial_monotonic.all()),
      "; each trial one contiguous run:", bool((vr.n_trials == vr.n_trial_runs).all()))
print("trial type codes present:", sorted({t for s in vr.trial_types_present for t in s}))
print("position range (cm): min", round(vr.x_min.min(), 4), "max", round(vr.x_max.max(), 2))
v200 = vr[vr.track_length == c.TRACK_LENGTH_CM]
print(f"position range (cm) in the {len(v200)} sessions on the 200 cm track: min", round(v200.x_min.min(), 4),
      "max", round(v200.x_max.max(), 2), "; per-session maximum from", round(v200.x_max.min(), 2))
print("\nTrack length (median of the per-trial maximum position, to 10 cm), sessions by mouse:")
print(pd.crosstab(vr.mouse, vr.track_length, margins=True).to_string())
gap = vr.t_last - vr.es_duration
print(f"\nPosition end minus raw signal end (s): median {gap.median():.4f}, "
      f"within 0.1 s in {(gap.abs() < 0.1).sum()} of {len(vr)} files")
print(vr.loc[gap.abs() >= 0.1, ["path", "es_duration", "t_last", "n_trials"]].round(1).to_string(index=False))
print("\nExample trial type sequences (first 30 trials):")
print(vr.groupby("mouse").trial_type_sequence.first().to_string())
print("\nSessions with probe trials:", int((vr.n_probe > 0).sum()), "of", len(vr))
print("median speed over each session (cm/s): range", round(vr.median_speed.min(), 1), "to", round(vr.median_speed.max(), 1))

VR position: sample interval (s) [np.float64(0.001)] ; first timestamp [np.float64(0.0)]
trial numbers start at 1 in every file: True ; never decrease: True ; each trial one contiguous run: True
trial type codes present: [0, 1, 2]
position range (cm): min 0.0 max 342.0
position range (cm) in the 326 sessions on the 200 cm track: min 0.0 max 200.0 ; per-session maximum from 199.99

Track length (median of the per-trial maximum position, to 10 cm), sessions by mouse:
track_length  200.0  230.0  270.0  340.0  All
mouse                                        
M1               11      0      0      0   11
M10              25      0      0      0   25
M11              36      2      5      2   45
M12              28      0      2      0   30
M13              29      0      0      0   29
M14              37      1      5      2   45
M15              21      0      0      0   21
M2               29      0      0      0   29
M3               30      0      5      0   35
M4               30      0      0      0   30
M6               22      0      0      0   22
M7               28      0      5      0   33
All             326      3     22      4  355

Position end minus raw signal end (s): median -0.0007, within 0.1 s in 353 of 355 files
                                                       path  es_duration  t_last  n_trials
sub-M13/sub-M13_ses-M13-D27-2021-06-15-11-43-42_ecephys.nwb       2167.6   687.5      21.0
 sub-M15/sub-M15_ses-M15-D5-2021-05-14-12-46-13_ecephys.nwb       1778.6   881.4      14.0

Example trial type sequences (first 30 trials):
mouse
M1     NBBBBNBBBBNBBBBNBBBBNBBBBNBBBB
M10       BNBBNBBNBBNBBNBBNBBNBBNBBNB
M11                           BBBNBBN
M12    NBBPBBNBBNBBPBBNBBNBBPBBNBBNBB
M13               BBBNNBBBBBNBBBNBBNB
M14            BNNBBBNBBBBNBBBBNBBNNN
M15                          BBBNBBNB
M2         BBNBBBBBNBBBBNBBBBNBBBBNBB
M3     NBBBBPBBBBBNBBBBPBBBBNBBBBBPBB
M4                            BBBBNBB
M6     BNBBBPBBBNBBBPBBBNBBBPBBBNBBBP
M7     BBBBNBBBBNBBBBNBBBBNBBBBBBNBBB

Sessions with probe trials: 105 of 355
median speed over each session (cm/s): range 0.0 to 70.7

# The paper's summary table (Supplementary file 1): sessions and trials per session per mouse.
paper_table = pd.DataFrame({
    "mouse": ["M1", "M3", "M6", "M7", "M10", "M11", "M12", "M13", "M14"],
    "paper_sessions": [11, 26, 21, 23, 16, 35, 21, 25, 37],
    "paper_trials_mean": [45.5, 270.3, 172.3, 268.1, 74.3, 173.1, 244.4, 47.5, 316.4],
    "paper_trials_sd": [33.6, 93.5, 99.8, 137.8, 34.8, 126.1, 68.5, 25.1, 148.2],
}).set_index("mouse")
# Sessions the paper's code can use: 200 cm track and at least 10 trials (claims file, opening comment).
use = vr[(vr.track_length == c.TRACK_LENGTH_CM) & (vr.n_trials >= c.MIN_TRIALS)]
mine = use.groupby("mouse").n_trials.agg(sessions="size", trials_mean="mean", trials_sd="std")
# The paper sorted open arena and VR recordings of a day together, so a VR session without an open arena
# session on the same day cannot give grid cells.
paired_days = set(of.groupby(["mouse", "day"]).size().index)
use_paired = use[[(m, d) in paired_days for m, d in zip(use.mouse, use.day)]]
mine["sessions_with_open_arena_same_day"] = use_paired.groupby("mouse").size()
cmp = paper_table.join(mine).round(1)
print(cmp.to_string())
print(f"\nNine mice: paper {paper_table.paper_sessions.sum()} sessions; dandiset {int(cmp.sessions.sum())} "
      f"200 cm sessions with >= 10 trials, {int(cmp.sessions_with_open_arena_same_day.sum())} of them with an open arena session the same day")
print("Mean over the nine mice of each mouse's mean trials per session: paper 179.1; dandiset",
      round(cmp.trials_mean.mean(), 1))
others = mine.drop(paper_table.index, errors="ignore")
print("\nMice in the dandiset but not in the paper's table (the code drops M2, M4 and M15, which had no grid cells):")
print(others.round(1).to_string())
print("\nDays by session types held (o = open arena, v = VR):")
print(survey.groupby(["mouse", "day"]).session_description.agg(lambda s: "".join(sorted(x[0] for x in s))).value_counts().to_string())

       paper_sessions  paper_trials_mean  paper_trials_sd  sessions  trials_mean  trials_sd  sessions_with_open_arena_same_day
mouse                                                                                                                         
M1                 11               45.5             33.6        11         58.0       38.2                                 11
M3                 26              270.3             93.5        30        259.2      109.1                                 29
M6                 21              172.3             99.8        21        183.3      102.1                                 21
M7                 23              268.1            137.8        28        271.7      144.3                                 27
M10                16               74.3             34.8        25        106.5       74.3                                 25
M11                35              173.1            126.1        35        165.1      127.7                                 35
M12                21              244.4             68.5        28        288.1       86.8                                 28
M13                25               47.5             25.1        27         44.5       24.1                                 26
M14                37              316.4            148.2        37        278.9      155.8                                 37

Nine mice: paper 215 sessions; dandiset 242 200 cm sessions with >= 10 trials, 239 of them with an open arena session the same day
Mean over the nine mice of each mouse's mean trials per session: paper 179.1; dandiset 183.9

Mice in the dandiset but not in the paper's table (the code drops M2, M4 and M15, which had no grid cells):
       sessions  trials_mean  trials_sd  sessions_with_open_arena_same_day
mouse                                                                     
M15          12         31.4       12.7                                 12
M2           29         78.9       46.7                                 21
M4           29        189.7       86.4                                 28

Days by session types held (o = open arena, v = VR):
session_description
ov    331
o      48
v      24

# Open arena sessions: nine undescribed columns. Their ranges suggest what they are.
cmin = pd.DataFrame(of.col_min.tolist()); cmax = pd.DataFrame(of.col_max.tolist())
print("Open arena position: sample interval (s)", sorted(of.dt.round(6).unique()))
print("Per-column range over sessions, median of the session minimum and maximum:")
print(pd.DataFrame({"min_min": cmin.min(), "min_median": cmin.median(), "max_median": cmax.median(), "max_min": cmax.min(),
                    "max_max": cmax.max()}).round(2).to_string())
print("NaN fraction, largest over sessions, per column:", np.round(np.array(of.col_nan.tolist()).max(axis=0), 3).tolist())
of["x_extent"] = cmax[1].values - cmin[1].values
of["y_extent"] = cmax[3].values - cmin[3].values
small = of[(of.x_extent < 90) & (of.y_extent < 90)]
print(f"\nSessions whose columns 1 and 3 span less than 90 cm in both: {len(small)}, by mouse:",
      small.mouse.value_counts().to_dict())
print("their extents (cm): x", round(small.x_extent.min(), 1), "to", round(small.x_extent.max(), 1),
      "; y", round(small.y_extent.min(), 1), "to", round(small.y_extent.max(), 1))
rest = of.drop(small.index)
print("other sessions: x", round(rest.x_extent.min(), 1), "to", round(rest.x_extent.max(), 1),
      "; y", round(rest.y_extent.min(), 1), "to", round(rest.y_extent.max(), 1))
gap = of.t_last - of.es_duration
print(f"\nPosition end minus raw signal end (s): median {gap.median():.2f}, 5th to 95th percentile "
      f"{gap.quantile(0.05):.2f} to {gap.quantile(0.95):.2f}; more than 10 s apart in {(gap.abs() > 10).sum()} of {len(of)} files"
      f" (largest {gap.max():.1f} s longer, {-gap.min():.1f} s shorter)")

Open arena position: sample interval (s) [np.float64(0.033333)]
Per-column range over sessions, median of the session minimum and maximum:
   min_min  min_median  max_median  max_min   max_max
0     0.00        1.13     1816.73  1409.80   3225.80
1     0.00        0.00      101.99    67.50    103.69
2     0.00        0.00      448.77   296.99    456.23
3     0.00        0.00      100.25    67.63    103.20
4     0.00        0.00      441.11   297.57    454.09
5  -180.00     -179.99      180.00   179.95    180.00
6    33.43     1793.02     5656.57  1983.43  37697.43
7     0.00        0.00       70.48    32.09    131.62
8   -17.93        0.07     1814.07  1402.80   3225.23
NaN fraction, largest over sessions, per column: [0.0, 0.008, 0.008, 0.008, 0.008, 0.008, 0.0, 0.008, 0.0]

Sessions whose columns 1 and 3 span less than 90 cm in both: 38, by mouse: {'M2': 20, 'M1': 18}
their extents (cm): x 67.5 to 89.4 ; y 67.6 to 89.8
other sessions: x 89.3 to 103.7 ; y 88.0 to 103.2

Position end minus raw signal end (s): median 2.87, 5th to 95th percentile 1.20 to 9.49; more than 10 s apart in 20 of 379 files (largest 451.9 s longer, 178.2 s shorter)

# Recordings deposited twice: identical position series, or identical raw blocks in files with the same
# number of raw samples.
print("VR files with an identical position series:", int(vr.md5.duplicated(keep=False).sum()))
print("open arena files with an identical position series:", int(of.md5.duplicated(keep=False).sum()))
same_len = survey[survey.es_samples.duplicated(keep=False)]
print(f"files sharing a raw sample count with another file: {len(same_len)} "
      f"({same_len.es_samples.nunique()} distinct counts)")
digests = pd.DataFrame(c.map_files(c.raw_block_digest, [a for a in assets if a["path"] in set(same_len.path)], workers=16))
same_len = same_len.merge(digests, on="path", suffixes=("", "_raw"))
dup_raw = same_len[same_len.duplicated(["es_samples", "md5_raw"], keep=False)]
print("of these, files whose 1 s raw block at sample 3,000,000 matches another file's:", len(dup_raw))
print("median SD of those raw blocks (ADC units):", round(float(same_len["std"].median()), 1))
print("same-count pairs within one mouse:", int(same_len.groupby("es_samples").mouse.nunique().eq(1).sum()))

VR files with an identical position series: 0
open arena files with an identical position series: 0
files sharing a raw sample count with another file: 74 (36 distinct counts)

of these, files whose 1 s raw block at sample 3,000,000 matches another file's: 0
median SD of those raw blocks (ADC units): 590.9
same-count pairs within one mouse: 2

rng = np.random.default_rng(943)
check_vr = []
for m in ["M3", "M7", "M11", "M12", "M13", "M14"]:
    pool = use[use.mouse == m].sort_values("path")
    check_vr.append(pool.iloc[rng.integers(len(pool))].path)
pool_of = of[~of.mouse.isin(["M1", "M2"])].sort_values("path")
check_of = list(pool_of.iloc[rng.choice(len(pool_of), 3, replace=False)].path)
t = time.time()
mua = c.map_files(c.mua_check, [by_path[p] for p in check_vr + check_of], workers=3)
print(f"read and filtered {len(mua)} five-minute raw blocks in {time.time() - t:.0f} s")
shifts = np.array([-50, -20, -10, -5, -2, -1, 0, 1, 2, 5, 10, 20, 50])
rows = []
for p, r in zip(check_vr + check_of, mua):
    lags, cc = c.lagged_corr(r["mua"], r["speed"], 50)
    con = c.shifted_contrast(r["mua"], r["speed"], shifts)
    rows.append({"file": p.split("/")[1][:-12], "type": "vr" if p in check_vr else "openfield",
                 "r_lag0": cc[50], "peak_lag_s": lags[np.argmax(cc)] * 0.1, "r_peak": cc.max(),
                 "r_at_+-5s": np.round([cc[0], cc[-1]], 3).tolist(),
                 "run/stop_lag0": con[shifts == 0][0], "best_shift_s": shifts[np.argmax(con)] * 0.1,
                 "run/stop_at_+-5s": np.round([con[0], con[-1]], 3).tolist(),
                 "stopped_bins": int((r["speed"] < c.STOP_SPEED_CM_S).sum()), "bins": len(r["speed"])})
clock = pd.DataFrame(rows)
print(clock.round(3).to_string(index=False))

read and filtered 9 five-minute raw blocks in 43 s
                                   file      type  r_lag0  peak_lag_s  r_peak       r_at_+-5s  run/stop_lag0  best_shift_s run/stop_at_+-5s  stopped_bins  bins
  sub-M3_ses-M3-D27-2020-12-04-15-10-38        vr  -0.198         1.8   0.302  [0.22, -0.174]          0.617           2.0   [1.291, 0.992]           536  3000
   sub-M7_ses-M7-D7-2020-11-06-16-14-44        vr   0.200         0.1   0.207  [0.111, 0.152]          1.192           0.1   [1.087, 1.135]          1578  3000
sub-M11_ses-M11-D12-2021-05-25-09-49-23        vr   0.123        -0.3   0.147 [-0.014, 0.036]          1.219          -0.1    [0.984, 1.04]          2331  3000
sub-M12_ses-M12-D17-2021-06-01-11-08-39        vr  -0.134         1.7   0.193 [-0.049, -0.06]          0.791           2.0    [0.876, 0.92]           738  3000
sub-M13_ses-M13-D13-2021-05-26-10-56-26        vr   0.093         0.0   0.093  [0.005, 0.014]          1.066          -0.1   [1.001, 1.002]          1410  3000
sub-M14_ses-M14-D34-2021-06-24-12-48-57        vr   0.129         0.4   0.149   [0.009, 0.04]          1.002           1.0   [1.028, 1.036]           537  3000
sub-M11_ses-M11-D35-2021-06-25-11-26-39 openfield   0.100         0.4   0.193 [-0.004, 0.033]          1.063           0.1   [1.002, 1.007]           617  3000
sub-M14_ses-M14-D22-2021-06-08-11-54-13 openfield   0.249         0.0   0.249  [0.036, 0.053]          1.173          -0.1   [1.065, 1.017]           485  3000
  sub-M3_ses-M3-D19-2020-11-22-15-25-38 openfield   0.327         0.3   0.343 [-0.005, 0.043]          1.272           0.0   [1.007, 1.003]          1003  3000

fig, axes = plt.subplots(1, 2, figsize=(11, 3.5))
for p, r in zip(check_vr + check_of, mua):
    lags, cc = c.lagged_corr(r["mua"], r["speed"], 50)
    axes[0].plot(lags * 0.1, cc, lw=1, label=p.split("/")[1][4:-12], ls="-" if p in check_vr else "--")
    axes[1].plot(shifts * 0.1, c.shifted_contrast(r["mua"], r["speed"], shifts), lw=1, marker=".")
axes[0].set(xlabel="lag of speed relative to multi-unit activity (s)", ylabel="Pearson r", title="MUA against speed")
axes[1].set(xlabel="shift of running/stopped labels (s)", ylabel="MUA running / stopped", title="Shifted labels")
axes[0].legend(fontsize=6)
plt.tight_layout(); plt.show()
print("plotted ranges: r from", round(min(c.lagged_corr(r['mua'], r['speed'], 50)[1].min() for r in mua), 3),
      "to", round(max(c.lagged_corr(r['mua'], r['speed'], 50)[1].max() for r in mua), 3))
print("sessions whose correlation peaks within 0.5 s of zero lag:", int((clock.peak_lag_s.abs() <= 0.5).sum()), "of", len(clock))
print("sessions with more activity running than stopped at no shift:", int((clock["run/stop_lag0"] > 1).sum()), "of", len(clock))
print("sessions whose running/stopped contrast is largest within 0.2 s of no shift:", int((clock.best_shift_s.abs() <= 0.2).sum()), "of", len(clock))

plotted ranges: r from -0.278 to 0.343
sessions whose correlation peaks within 0.5 s of zero lag: 7 of 9
sessions with more activity running than stopped at no shift: 7 of 9
sessions whose running/stopped contrast is largest within 0.2 s of no shift: 6 of 9

rate = c.streaming_rate([by_path[p] for p in use.sort_values("path").path.iloc[::20][:16]], workers=16)
print(f"read {rate['bytes'] / 1e6:.0f} MB in {rate['wall_s']:.1f} s: aggregate {rate['aggregate_mb_s']:.1f} MB/s "
      f"with {rate['workers']} processes; median single stream {rate['single_stream_mb_s']:.1f} MB/s")
nine = ["M1", "M3", "M6", "M7", "M10", "M11", "M12", "M13", "M14"]
days = use_paired[use_paired.mouse.isin(nine)][["mouse", "day"]].drop_duplicates()
need = survey.merge(days, on=["mouse", "day"])
gb_all = survey.size_gb.sum()
gb_need = need.size_gb.sum()
for label, gb in [("every file", gb_all), ("days of the nine mice with a usable VR session and an open arena session", gb_need)]:
    hours = gb * 1e3 / rate["aggregate_mb_s"] / 3600
    print(f"{label}: {gb:.0f} GB, {len(survey) if gb == gb_all else len(need)} files; reading alone at the measured rate "
          f"takes {hours:.1f} h (notebook limit {c.NOTEBOOK_LIMIT_MIN} min)")
# These days hold the 239 paired sessions, a superset of the paper's 215; scaled to 215 as a smaller bound.
n_paired = int(use_paired.mouse.isin(nine).sum())
gb_scaled = gb_need * paper_table.paper_sessions.sum() / n_paired
print(f"scaled by the paper's {paper_table.paper_sessions.sum()} of these {n_paired} sessions: {gb_scaled:.0f} GB; "
      f"reading alone takes {gb_scaled * 1e3 / rate['aggregate_mb_s'] / 3600:.1f} h")
print(f"\nnotebook so far: {(time.time() - started) / 60:.1f} min")

read 307 MB in 3.2 s: aggregate 95.8 MB/s with 16 processes; median single stream 9.8 MB/s
every file: 976 GB, 734 files; reading alone at the measured rate takes 2.8 h (notebook limit 60 min)
days of the nine mice with a usable VR session and an open arena session: 638 GB, 478 files; reading alone at the measured rate takes 1.8 h (notebook limit 60 min)
scaled by the paper's 215 of these 239 sessions: 574 GB; reading alone takes 1.7 h

notebook so far: 2.3 min

Claim	Status
grid-cell-count, supp1-recording-summary, grid-field-size-spacing	blocked, too_large (needs sorted units)
fig2-session-level-classes, fig2e-nongrid-classes, fig2f-spatial-information, fig2f-firing-rate-similar, fig2f-periodogram-features-similar, fig2f-aperiodic-low-rate	blocked, too_large
fig3f-switching-within-session, fig3f-session-class-bias, fig3g-trial-spatial-information, fig3g-trial-firing-rate-similar, fig3s4-stable-cells, fig3s5-grid-vs-nongrid-modes, fig3s6-block-lengths, fig3s7-trial-type-ratio	blocked, too_large
fig4e-grid-agreement-higher, fig4e-agreement-above-shuffle, fig4f-nongrid-spatial-information	blocked, too_large
fig6b-stable-beaconed-similar, fig6b-stable-uncued-reduced, fig6a-stable-stop-density, fig6d-epochs-beaconed-similar, fig6d-epochs-nonbeaconed-reduced, fig6d-epochs-probe, fig6c-epoch-stop-density, fig6e-hit-try-run, fig6f-odds-ratio, fig6s1-speed-profiles, fig6s2-template-method	blocked, too_large

Clark and Nolan (2024) on dandiset 000943¶

Data check¶

Shared clock of the raw signal and the position series¶

Why every claim is blocked: the size of spike sorting¶

Summary¶