Reproducing Guan et al. (preprint 2022) from dandiset 000252¶
Paper. Guan, Aflalo, Kadlec, Gámez de Leon, Rosario, Bari, Pouratian and Andersen, "Compositional coding of individual finger movements in human posterior parietal cortex and motor cortex enables ten-finger decoding", medRxiv/bioRxiv preprint, doi:10.1101/2022.12.07.22283227.
Dataset. DANDI:000252, "Finger_RL: human intracortical recordings during attempted finger movements of right and left hands", version 0.230408.2207 (doi:10.48324/dandi.000252/0.230408.2207). The paper cites no version (its Data Availability gives only the dandiset URL), and 0.230408.2207 is the only published version, so it is pinned. The draft was modified later (2023-06-20) but lists the same 12 assets with the same asset identifiers; the check below prints this.
Text the claims file was written from. The bioRxiv preprint, version 1, full text with Methods and the supplementary figure legends (claims file claims/000252/10.1101-2022.12.07.22283227.yaml). The paper links no analysis code.
What the dataset holds. Only the ten-finger press task: 10 sessions of participant N (subject P1 in the files) and 2 of participant J (subject N1). The alternating-cues task (Figures 4 and 5) and the BMI tasks (Figure 6, Supplementary Figures 6 and 7) are not in it.
| Claim | Tested here |
|---|---|
| tenfinger-recording-counts | tested |
| tenfinger-pooled-npcc-neurons | tested |
| tenfinger-discriminative-fraction | tested |
| tenfinger-decoding-accuracy | tested |
| tenfinger-error-structure | tested |
| fig8b-matching-closer | tested |
| fig8a-npcc-index-distinct | tested |
| fig8e-hand-ccgp | tested |
| fig8e-finger-ccgp | tested |
| jppc-finger-ccgp-low | tested |
| figs11-jppc-lower-matched-size | tested |
| fig4c-contra-discriminative | blocked, data_missing (alternating-cues task not in the dataset) |
| fig4e-plan-vs-move | blocked, data_missing |
| fig5a-contra-decoding | blocked, data_missing |
| fig5bcd-npcc-data-requirements | blocked, data_missing |
| fig5e-planning-decoding | blocked, data_missing |
| fig6-online-bmi | blocked, data_missing (BMI sessions not in the dataset) |
| figs6-single-array-retrospective | blocked, data_missing |
| figs7-left-hand-bmi | blocked, data_missing |
The readings O1 to O11 named below are those of the claims file's header.
import os, time, warnings
os.environ.setdefault("OMP_NUM_THREADS", "1")
warnings.filterwarnings("ignore")
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy.stats import bootstrap
import requests
import fingers as F
pd.set_option("display.width", 200)
BLUE, ORANGE = "#2a78d6", "#eb6834"
T0 = time.time()
print("dandiset", F.DANDISET_ID, "version", F.VERSION)
dandiset 000252 version 0.230408.2207
Data check¶
meta = F.version_metadata()
print("title:", meta["name"])
print("related:", [(r.get("relation"), r.get("identifier")) for r in meta.get("relatedResource", [])])
assets = F.list_assets()
draft = requests.get(f"{F.API}/dandisets/{F.DANDISET_ID}/versions/draft/assets/?page_size=200").json()["results"]
print("assets in version:", len(assets), "| in draft:", len(draft),
"| same asset ids:", {a["asset_id"] for a in assets} == {a["asset_id"] for a in draft})
t = time.time()
S = [F.load_session(a) for a in assets]
print(f"loaded {len(S)} files in {time.time() - t:.0f} s")
rows = []
for s in S:
tr = s.trials
rows.append(dict(
file=s.path.split("/")[-1], subject=s.subject, sex=s.sex, date=s.date, task=s.description,
trials=len(tr), runs=tr.Run.nunique(), per_condition=f"{tr.finger.value_counts().min()}-{tr.finger.value_counts().max()}",
units_PPC=int((s.group == "PPC").sum()), units_MC=int((s.group == "MC").sum()),
electrodes=s.extra["electrode_groups"], unsorted_units=s.n_unsorted, spikes_outside_obs=s.n_outside_obs))
summary = pd.DataFrame(rows)
print(summary.to_string())
for s in S[:1] + S[2:3]:
print(s.subject, s.sex, s.extra["electrode_locations"])
title: Finger_RL: human intracortical recordings during attempted finger movements of right and left hands
related: [('dcite:IsDescribedBy', 'https://doi.org/10.1101/2022.12.07.22283227 '), ('dcite:Continues', 'DANDI:000147')]
assets in version: 12 | in draft: 12 | same asset ids: True
loaded 12 files in 10 s
file subject sex date task trials runs per_condition units_PPC units_MC electrodes unsorted_units spikes_outside_obs
0 sub-N1_ses-2019-02-12_ecephys.nwb N1 M 2019-02-12 ten-finger press task, with delay 100 2 10-10 159 139 {'PPC': 96, 'MC': 95} 283 0
1 sub-N1_ses-2019-12-20_ecephys.nwb N1 M 2019-12-20 ten-finger press task, with delay 100 2 10-10 160 120 {'MC': 96, 'PPC': 96} 224 0
2 sub-P1_ses-2015-08-12_ecephys.nwb P1 F 2015-08-12 ten-finger press task, with delay 100 2 10-10 112 0 {'PPC': 96} 101 0
3 sub-P1_ses-2015-08-14_ecephys.nwb P1 F 2015-08-14 ten-finger press task, with delay 100 2 10-10 102 0 {'PPC': 96} 92 0
4 sub-P1_ses-2015-08-17_ecephys.nwb P1 F 2015-08-17 ten-finger press task, with delay 100 2 10-10 119 0 {'PPC': 96} 93 0
5 sub-P1_ses-2015-08-21_ecephys.nwb P1 F 2015-08-21 ten-finger press task, with delay 100 2 10-10 110 0 {'PPC': 96} 94 0
6 sub-P1_ses-2016-04-27_ecephys.nwb P1 F 2016-04-27 ten-finger press task, with delay 100 2 10-10 116 0 {'PPC': 96} 108 0
7 sub-P1_ses-2016-05-04_ecephys.nwb P1 F 2016-05-04 ten-finger press task, with delay 100 2 10-10 110 0 {'PPC': 96} 108 0
8 sub-P1_ses-2016-05-11_ecephys.nwb P1 F 2016-05-11 ten-finger press task, with delay 100 2 10-10 118 0 {'PPC': 96} 113 0
9 sub-P1_ses-2016-05-16_ecephys.nwb P1 F 2016-05-16 ten-finger press task, with delay 100 2 10-10 109 0 {'PPC': 96} 96 0
10 sub-P1_ses-2016-05-20_ecephys.nwb P1 F 2016-05-20 ten-finger press task, with delay 100 2 10-10 105 0 {'PPC': 96} 97 0
11 sub-P1_ses-2016-05-23_ecephys.nwb P1 F 2016-05-23 ten-finger press task, with delay 100 2 10-10 113 0 {'PPC': 96} 100 0
N1 M ['left motor cortex (MC), near hand knob', 'left posterior parietal cortex (PPC), superior parietal lobule (SPL)']
P1 F ['left posterior parietal cortex (PPC), junction of postcentral and intraparietal sulci (PC-IP)']
# Conditions, phase durations and units of the variables used
s0 = S[0]
print("conditions:", sorted(s0.trials.finger.unique()), "| hand:", sorted(s0.trials.hand.unique()),
"| finger_type:", sorted(s0.trials.finger_type.unique()))
print("finger == hand + finger_type in every trial:",
all((s.trials.finger == s.trials.hand + s.trials.finger_type).all() for s in S))
dur = []
for s in S:
tr = s.trials
dur.append(dict(file=s.path.split("/")[-1],
iti=(tr.cue_on_time - tr.start_time).median(), cue=(tr.cue_off_time - tr.cue_on_time).median(),
delay_min=(tr.go_on_time - tr.cue_off_time).min(), delay_max=(tr.go_on_time - tr.cue_off_time).max(),
go=(tr.stop_time - tr.go_on_time).median(), go_min=(tr.stop_time - tr.go_on_time).min(),
ordered=bool(((tr.start_time <= tr.cue_on_time) & (tr.cue_on_time < tr.cue_off_time) &
(tr.cue_off_time < tr.go_on_time) & (tr.go_on_time < tr.stop_time)).all())))
print("phase durations, s (median unless stated):")
print(pd.DataFrame(dur).round(2).to_string())
allrates = np.concatenate([F.rates(s, "go").ravel() for s in S])
print(f"Go-window rates (Hz): min {allrates.min():.1f}, median {np.median(allrates):.1f}, max {allrates.max():.1f}")
print("spike-time range per file (s):", [(round(min(x.min() for x in s.spikes if len(x))), round(max(x.max() for x in s.spikes if len(x)))) for s in S])
conditions: ['Li', 'Lm', 'Lp', 'Lr', 'Lt', 'Ri', 'Rm', 'Rp', 'Rr', 'Rt'] | hand: ['L', 'R'] | finger_type: ['i', 'm', 'p', 'r', 't']
finger == hand + finger_type in every trial: True
phase durations, s (median unless stated):
file iti cue delay_min delay_max go go_min ordered
0 sub-N1_ses-2019-02-12_ecephys.nwb 1.57 1.07 1.06 1.58 1.07 1.07 True
1 sub-N1_ses-2019-12-20_ecephys.nwb 2.09 1.61 0.86 1.14 2.08 2.05 True
2 sub-P1_ses-2015-08-12_ecephys.nwb 2.10 1.61 1.08 1.15 1.12 1.09 True
3 sub-P1_ses-2015-08-14_ecephys.nwb 2.10 1.61 1.08 1.13 1.11 1.09 True
4 sub-P1_ses-2015-08-17_ecephys.nwb 2.10 1.60 1.08 1.12 1.12 1.10 True
5 sub-P1_ses-2015-08-21_ecephys.nwb 2.10 1.60 1.08 1.13 1.11 1.09 True
6 sub-P1_ses-2016-04-27_ecephys.nwb 1.59 1.10 1.08 1.62 1.12 1.09 True
7 sub-P1_ses-2016-05-04_ecephys.nwb 1.59 1.10 1.09 1.61 1.12 1.09 True
8 sub-P1_ses-2016-05-11_ecephys.nwb 1.59 1.11 1.07 1.59 1.12 1.09 True
9 sub-P1_ses-2016-05-16_ecephys.nwb 1.59 1.10 1.11 1.61 1.12 1.09 True
10 sub-P1_ses-2016-05-20_ecephys.nwb 1.59 1.11 1.09 1.61 1.12 1.09 True
11 sub-P1_ses-2016-05-23_ecephys.nwb 1.59 1.11 1.09 1.62 1.12 1.09 True
Go-window rates (Hz): min 0.0, median 2.0, max 118.0
spike-time range per file (s): [(171, 736), (156, 1031), (276, 1021), (310, 1054), (132, 814), (168, 819), (246, 1025), (26, 577), (65, 859), (49, 606), (21, 726), (23, 629)]
The files name participant N P1 (female, one PPC array at PC-IP) and participant J N1 (male, PPC array in SPL and MC array near the hand knob), as the sex and electrode locations above show. The subject identifier N1 is participant J, not N. Every Go phase lasts at least the 0.7 s the Go window needs.
# Recordings deposited twice: compare trial times and spike trains between all pairs of files
dups = []
for i in range(len(S)):
for j in range(i + 1, len(S)):
a, b = S[i], S[j]
same_trials = len(a.trials) == len(b.trials) and np.allclose(a.trials.go_on_time.values, b.trials.go_on_time.values)
sa = {tuple(np.round(x[:30], 4)) for x in a.spikes}
sb = {tuple(np.round(x[:30], 4)) for x in b.spikes}
if same_trials or (sa & sb):
dups.append((a.path, b.path, same_trials, len(sa & sb)))
print("pairs of files with identical trial times or identical spike trains:", dups if dups else "none")
within = sum(len(s.spikes) - len({tuple(np.round(x[:30], 4)) for x in s.spikes}) for s in S)
print("units with an identical spike train to another unit of the same file:", within)
pairs of files with identical trial times or identical spike trains: none units with an identical spike train to another unit of the same file: 0
# Clock check: decode the cued finger in 500 ms windows starting from -3 s to +1.5 s around the Go cue,
# The ITI should be at chance; decoding should rise after the cue and Go. The shuffled-label control is per session, below.
rng = np.random.default_rng(1)
starts = np.round(np.arange(-3.0, 1.51, 0.5), 2)
tc = {}
for arr in F.ARRAYS:
accs, sh = [], []
for t0 in starts:
T, P = [], []
for s in [s for s in S if s.subject == F.ARRAYS[arr][0]]:
m = (s.group == F.ARRAYS[arr][1]) & F.keep_mask(s, "a")
X = F.counts(s.spikes, s.trials.go_on_time.values, t0, 0.5)[:, m] / 0.5
y = s.trials.finger.values
T.append(y); P.append(F.cv_predict(X, y))
T = np.concatenate(T)
accs.append(np.mean(T == np.concatenate(P)))
tc[arr] = accs
print("window start rel. Go (s):", list(starts))
for arr, a in tc.items():
print(f"{arr:6s} accuracy:", np.round(a, 2).tolist())
# per-session ITI and Go accuracy, and with labels shuffled
for s in S:
for arr in [a for a in F.ARRAYS if F.ARRAYS[a][0] == s.subject]:
X_go, y = F.array_data(s, arr, "a", "go"); X_iti, _ = F.array_data(s, arr, "a", "iti")
print(f"{s.subject} {s.date} {arr:6s} Go {np.mean(F.cv_predict(X_go, y) == y):.2f} ITI {np.mean(F.cv_predict(X_iti, y) == y):.2f}"
f" Go with shuffled labels {np.mean(F.cv_predict(X_go, rng.permutation(y)) == y):.2f} (chance 0.10)")
window start rel. Go (s): [np.float64(-3.0), np.float64(-2.5), np.float64(-2.0), np.float64(-1.5), np.float64(-1.0), np.float64(-0.5), np.float64(0.0), np.float64(0.5), np.float64(1.0), np.float64(1.5)] N-PPC accuracy: [0.1, 0.12, 0.27, 0.37, 0.41, 0.51, 0.66, 0.68, 0.47, 0.33] J-PPC accuracy: [0.04, 0.15, 0.23, 0.4, 0.56, 0.68, 0.74, 0.61, 0.52, 0.42] J-MC accuracy: [0.06, 0.12, 0.26, 0.44, 0.62, 0.71, 0.82, 0.82, 0.66, 0.52] N1 2019-02-12 J-PPC Go 0.76 ITI 0.05 Go with shuffled labels 0.09 (chance 0.10) N1 2019-02-12 J-MC Go 0.79 ITI 0.03 Go with shuffled labels 0.06 (chance 0.10) N1 2019-12-20 J-PPC Go 0.63 ITI 0.12 Go with shuffled labels 0.11 (chance 0.10) N1 2019-12-20 J-MC Go 0.78 ITI 0.07 Go with shuffled labels 0.08 (chance 0.10) P1 2015-08-12 N-PPC Go 0.72 ITI 0.08 Go with shuffled labels 0.06 (chance 0.10) P1 2015-08-14 N-PPC Go 0.66 ITI 0.17 Go with shuffled labels 0.04 (chance 0.10) P1 2015-08-17 N-PPC Go 0.70 ITI 0.08 Go with shuffled labels 0.11 (chance 0.10)
P1 2015-08-21 N-PPC Go 0.72 ITI 0.04 Go with shuffled labels 0.10 (chance 0.10) P1 2016-04-27 N-PPC Go 0.74 ITI 0.10 Go with shuffled labels 0.07 (chance 0.10) P1 2016-05-04 N-PPC Go 0.68 ITI 0.06 Go with shuffled labels 0.11 (chance 0.10) P1 2016-05-11 N-PPC Go 0.74 ITI 0.10 Go with shuffled labels 0.06 (chance 0.10) P1 2016-05-16 N-PPC Go 0.72 ITI 0.14 Go with shuffled labels 0.09 (chance 0.10) P1 2016-05-20 N-PPC Go 0.87 ITI 0.12 Go with shuffled labels 0.12 (chance 0.10) P1 2016-05-23 N-PPC Go 0.67 ITI 0.10 Go with shuffled labels 0.17 (chance 0.10)
Checks of the description against the data are stated after the outputs above: the condition labels, the per-session trial counts and the arrays are as the paper describes; the trial phases are in the order ITI, Cue, Delay, Go in every trial; Go-window decoding is high in every session while the ITI and shuffled labels are near the 10% chance, so spikes and trial table share a clock in each of the 12 sessions; the units' spike times are stored out of time order (sorted here before counting).
tenfinger-recording-counts and tenfinger-pooled-npcc-neurons¶
In the ten-finger task, N-PPC yielded 111 neurons per session (102 to 119) over 10 sessions of 100 trials, J-PPC 160 (159 to 160) and J-MC 130 (120 to 130) over 2 sessions of 100 trials.
1111 N-PPC neurons were recorded across the 10 sessions of the ten-finger task.
Counted before the 0.5 Hz exclusion and under both readings of O1: (a) mean rate over the unit's observed time (obs_intervals), (b) mean Go-window rate over all trials.
def sess_of(arr):
return [s for s in S if s.subject == F.ARRAYS[arr][0]]
count_rows = []
for o1 in ["none", "a", "b"]:
for arr in F.ARRAYS:
n = [int(F.keep_mask(s, o1)[s.group == F.ARRAYS[arr][1]].sum()) for s in sess_of(arr)]
count_rows.append(dict(O1=o1, array=arr, sessions=len(n), per_session=n, mean=np.mean(n), min=min(n), max=max(n), pooled=sum(n)))
counts_df = pd.DataFrame(count_rows)
print(counts_df.to_string())
print("trials per session:", sorted({len(s.trials) for s in S}))
O1 array sessions per_session mean min max pooled 0 none N-PPC 10 [112, 102, 119, 110, 116, 110, 118, 109, 105, 113] 111.4 102 119 1114 1 none J-PPC 2 [159, 160] 159.5 159 160 319 2 none J-MC 2 [139, 120] 129.5 120 139 259 3 a N-PPC 10 [95, 86, 94, 89, 100, 104, 100, 92, 92, 91] 94.3 86 104 943 4 a J-PPC 2 [157, 158] 157.5 157 158 315 5 a J-MC 2 [137, 115] 126.0 115 137 252 6 b N-PPC 10 [96, 91, 93, 82, 99, 100, 96, 91, 93, 89] 93.0 82 100 930 7 b J-PPC 2 [157, 155] 156.0 155 157 312 8 b J-MC 2 [137, 113] 125.0 113 137 250 trials per session: [100]
Result: reproduced, before the exclusion. The dataset has 10 N sessions and 2 J sessions of 100 trials each. Before the 0.5 Hz exclusion, N-PPC has a mean of 111.4 units per session (102 to 119), J-PPC 159.5 (159 to 160) and J-MC 129.5 (120 to 139). The J-MC maximum is 139 against the paper's 130, inside the 10% widening of that bound (143). After the exclusion the N-PPC mean falls to 94.3 (O1 a) or 93.0 (O1 b), outside tolerance, so the paper's counts are those before exclusion. Pooled over sessions, N-PPC has 1114 units before the exclusion (paper 1111), and 943 or 930 after it. The paper's counts thus match the units before the exclusion, while its Methods (Feature Extraction) apply the exclusion to all analyses; every analysis below therefore uses the units kept under O1 (a) or (b), as the claims file's readings say, and the counts before exclusion are only counts.
tenfinger-discriminative-fraction¶
Most neurons discriminate firing rates across the ten fingers (N-PPC 66%, J-PPC 57%, J-MC 78%, q < 0.05).
One-way ANOVA of Go-window rates across the ten conditions, Benjamini-Hochberg over all neurons of the array pooled across sessions (Methods, Single-neuron selectivity). Readings of O1. Variants, not readings of the criterion: no FDR correction, FDR per session, other analysis windows, and a sweep of the 500 ms window's start from 0.5 s before to 0.5 s after the Go cue (the Go phase lasts at least 1.07 s, so every window of the sweep ends before the trial does), with a 1 s window over the whole Go phase.
sel_rows = []
for o1 in ["none", "a", "b"]:
for arr in F.ARRAYS:
ps, per_sess = [], []
for s in sess_of(arr):
X, y = F.array_data(s, arr, o1)
p, _ = F.anova_eta2(X, y)
ps.append(p)
per_sess.append(np.mean(F.fdr_bh(p) < 0.05))
p = np.concatenate(ps)
q = F.fdr_bh(p)
sel_rows.append(dict(O1=o1, array=arr, neurons=len(p), constant_units=int(np.isnan(p).sum()),
pct_q05=100 * np.mean(q < 0.05), pct_uncorrected_p05=100 * np.mean(np.nan_to_num(p, nan=1) < 0.05),
pct_FDR_per_session_mean=100 * np.mean(per_sess),
per_session_pct=np.round(100 * np.array(per_sess)).astype(int).tolist()))
sel_df = pd.DataFrame(sel_rows)
print(sel_df.round(1).to_string())
print("paper: N-PPC 66, J-PPC 57, J-MC 78 (tolerance 10%)")
# variant: other windows (O1 a)
for arr in F.ARRAYS:
out = {}
for w in ["iti", "cue", "delay", "go"]:
p = np.concatenate([F.anova_eta2(*F.array_data(s, arr, "a", w))[0] for s in sess_of(arr)])
out[w] = round(100 * np.mean(F.fdr_bh(p) < 0.05), 1)
print(arr, "percent q < 0.05 by window (O1 a):", out)# variant: sweep of the window start around Go (O1 a), and a 1 s window from the Go cue
sweep = {}
for arr in F.ARRAYS:
out = {}
for st in np.round(np.arange(-0.5, 0.51, 0.1), 2):
p = np.concatenate([F.anova_eta2(*F.array_data(s, arr, "a", "go", start=st))[0] for s in sess_of(arr)])
out[f"{st + 0.0:+.1f}"] = round(100 * np.mean(F.fdr_bh(p) < 0.05), 1)
p = np.concatenate([F.anova_eta2(*F.array_data(s, arr, "a", "go", start=0.0, dur=1.0))[0] for s in sess_of(arr)])
out["0 to 1 s"] = round(100 * np.mean(F.fdr_bh(p) < 0.05), 1)
sweep[arr] = out
print("percent q < 0.05 by window start relative to Go (500 ms window; O1 a):")
print(pd.DataFrame(sweep).T.to_string())
O1 array neurons constant_units pct_q05 pct_uncorrected_p05 pct_FDR_per_session_mean per_session_pct
0 none N-PPC 1114 5 40.4 49.4 40.7 [41, 35, 45, 36, 38, 44, 38, 40, 55, 35]
1 none J-PPC 319 0 37.6 46.1 37.9 [39, 37]
2 none J-MC 259 2 63.7 66.8 62.3 [71, 53]
3 a N-PPC 943 0 47.8 55.0 47.4 [49, 42, 53, 44, 46, 46, 43, 47, 62, 42]
4 a J-PPC 315 0 38.1 46.3 38.4 [39, 37]
5 a J-MC 252 0 65.1 68.3 63.6 [72, 56]
6 b N-PPC 930 0 48.9 56.1 48.5 [49, 38, 54, 46, 46, 47, 49, 48, 61, 45]
7 b J-PPC 312 0 38.1 46.5 38.1 [38, 38]
8 b J-MC 250 0 65.6 68.8 64.1 [72, 57]
paper: N-PPC 66, J-PPC 57, J-MC 78 (tolerance 10%)
N-PPC percent q < 0.05 by window (O1 a): {'iti': np.float64(0.0), 'cue': np.float64(0.0), 'delay': np.float64(16.2), 'go': np.float64(47.8)}
J-PPC percent q < 0.05 by window (O1 a): {'iti': np.float64(0.0), 'cue': np.float64(0.0), 'delay': np.float64(36.8), 'go': np.float64(38.1)}
J-MC percent q < 0.05 by window (O1 a): {'iti': np.float64(0.0), 'cue': np.float64(0.8), 'delay': np.float64(49.2), 'go': np.float64(65.1)}
percent q < 0.05 by window start relative to Go (500 ms window; O1 a):
-0.5 -0.4 -0.3 -0.2 -0.1 +0.0 +0.1 +0.2 +0.3 +0.4 +0.5 0 to 1 s
N-PPC 30.6 34.7 36.3 39.2 40.1 41.1 45.6 47.8 45.6 42.1 42.7 60.8
J-PPC 43.5 41.9 43.5 41.9 42.9 40.6 42.2 38.1 32.1 34.0 29.2 56.8
J-MC 64.7 63.5 64.7 64.7 69.8 72.2 67.1 65.1 64.7 64.3 60.7 82.9
Result: not reproduced. Before the exclusion and under both readings of O1, the Go-window percentages of discriminative neurons are 40.4 to 48.9% for N-PPC, 37.6 to 38.1% for J-PPC and 63.7 to 65.6% for J-MC, against 66%, 57% and 78%; all three are below the lower tolerance bounds (59.4%, 51.3% and 70.2%). Without FDR correction they are still below the paper's values (at most 56.1%, 46.5% and 68.8%), FDR within each session, averaged over sessions, gives at most 0.3 points more than pooled FDR (for example 38.4% against 38.1% for J-PPC under O1 a), and none of the Delay, Cue and ITI windows gives more than the Go window. Moving the 500 ms window's start through the Go phase (from 0.5 s before to 0.5 s after the cue) gives at most 47.8% (N-PPC, start +0.2 s), 43.5% (J-PPC) and 72.2% (J-MC, start at the cue), so no 500 ms window brings N-PPC and J-PPC within tolerance. A 1 s window from the Go cue gives 60.8%, 56.8% and 82.9%, all within tolerance. The Methods (Feature Extraction) use 500 ms windows "except when otherwise specified", and the claims file fixes the 200 to 700 ms Go window (O11), so this is a variant only and does not change the outcome. It shows that the gap can be closed by a longer window, which the paper does not describe for this analysis. The ordering of the arrays (J-MC highest) is as in the paper.
tenfinger-decoding-accuracy and tenfinger-error-structure¶
A linear classifier discriminates all ten fingers well above the 10% chance level (cross-validated accuracy N-PPC 70%, J-PPC 66%, J-MC 75%).
Most ten-finger classification errors confuse an adjacent finger of the same hand or the matching finger of the other hand (N-PPC 76%, J-PPC 66%, J-MC 68% of errors).
Diagonal LDA per session, stratified 8-fold cross-validation, held-out predictions pooled over sessions (Methods, Offline classification). Readings O1 x O2; O2 (b) is averaged over seeds 0 to 9.
dec_rows, confusions = [], {}
for o1 in ["a", "b", "none"]:
for o2 in ["a", "b"]:
for arr in F.ARRAYS:
acc, share, per_sess = [], [], []
for seed in (range(10) if o2 == "b" else [0]):
T, P = [], []
for s in sess_of(arr):
X, y = F.array_data(s, arr, o1)
p = F.cv_predict(X, y, o2, seed)
T.append(y); P.append(p)
if seed == 0:
per_sess.append(round(np.mean(p == y), 2))
T, P = np.concatenate(T), np.concatenate(P)
e, adj, mat = F.error_kinds(T, P)
acc.append(np.mean(T == P)); share.append((adj + mat) / e)
if (o1, o2, seed) == ("a", "a", 0):
confusions[arr] = pd.crosstab(pd.Series(T, name="true"), pd.Series(P, name="pred")).reindex(index=F.CONDITIONS, columns=F.CONDITIONS, fill_value=0)
dec_rows.append(dict(O1=o1, O2=o2, array=arr, trials=len(T), accuracy=100 * np.mean(acc), acc_min_seed=100 * min(acc),
acc_max_seed=100 * max(acc), errors=e, adjacent=adj, matching=mat,
share_adj_or_match=100 * np.mean(share), per_session_acc=per_sess))
dec_df = pd.DataFrame(dec_rows)
print(dec_df.round(1).to_string())
print("paper accuracy: N-PPC 70, J-PPC 66, J-MC 75 (chance 10); error share: 76, 66, 68")
O1 O2 array trials accuracy acc_min_seed acc_max_seed errors adjacent matching share_adj_or_match per_session_acc 0 a a N-PPC 1000 72.2 72.2 72.2 278 119 94 76.6 [0.72, 0.66, 0.7, 0.72, 0.74, 0.68, 0.74, 0.72, 0.87, 0.67] 1 a a J-PPC 200 69.5 69.5 69.5 61 40 4 72.1 [0.76, 0.63] 2 a a J-MC 200 78.5 78.5 78.5 43 11 16 62.8 [0.79, 0.78] 3 a b N-PPC 1000 71.4 70.6 72.6 293 115 106 74.9 [0.77, 0.6, 0.71, 0.67, 0.71, 0.63, 0.77, 0.78, 0.83, 0.67] 4 a b J-PPC 200 68.8 66.5 72.0 56 35 2 67.1 [0.73, 0.64] 5 a b J-MC 200 78.4 77.0 80.0 42 13 14 64.6 [0.8, 0.79] 6 b a N-PPC 1000 73.0 73.0 73.0 270 108 98 76.3 [0.73, 0.69, 0.67, 0.77, 0.78, 0.63, 0.78, 0.74, 0.85, 0.66] 7 b a J-PPC 200 68.0 68.0 68.0 64 42 4 71.9 [0.74, 0.62] 8 b a J-MC 200 78.0 78.0 78.0 44 11 17 63.6 [0.77, 0.79] 9 b b N-PPC 1000 72.4 70.3 73.4 297 120 105 75.0 [0.78, 0.62, 0.72, 0.72, 0.75, 0.62, 0.8, 0.76, 0.82, 0.62] 10 b b J-PPC 200 68.2 66.0 72.0 56 34 2 66.2 [0.72, 0.62] 11 b b J-MC 200 78.0 76.5 79.5 43 13 15 65.1 [0.8, 0.79] 12 none a N-PPC 1000 71.7 71.7 71.7 283 114 98 74.9 [0.7, 0.64, 0.72, 0.7, 0.77, 0.66, 0.73, 0.73, 0.84, 0.68] 13 none a J-PPC 200 70.0 70.0 70.0 60 38 4 70.0 [0.75, 0.65] 14 none a J-MC 200 78.0 78.0 78.0 44 12 16 63.6 [0.78, 0.78] 15 none b N-PPC 1000 71.9 69.9 73.6 301 121 98 73.9 [0.78, 0.64, 0.73, 0.72, 0.71, 0.65, 0.76, 0.72, 0.84, 0.66] 16 none b J-PPC 200 69.2 67.5 72.5 55 34 2 67.6 [0.71, 0.65] 17 none b J-MC 200 78.4 77.0 80.0 42 14 14 65.3 [0.8, 0.79] paper accuracy: N-PPC 70, J-PPC 66, J-MC 75 (chance 10); error share: 76, 66, 68
fig, axes = plt.subplots(1, 3, figsize=(13, 4))
for ax, (arr, cm) in zip(axes, confusions.items()):
frac = cm.values / cm.values.sum(axis=1, keepdims=True)
im = ax.imshow(frac, cmap="Blues", vmin=0, vmax=1)
ax.set_xticks(range(10), F.CONDITIONS, fontsize=8); ax.set_yticks(range(10), F.CONDITIONS, fontsize=8)
ax.set_title(f"{arr}, O1 a, O2 a"); ax.set_xlabel("predicted"); ax.set_ylabel("cued")
fig.colorbar(im, ax=axes, shrink=0.8, label="fraction of cued trials")
plt.show()
print("confusion fractions are in [0, 1]; colour scale 0 to 1, no value outside it")
confusion fractions are in [0, 1]; colour scale 0 to 1, no value outside it
Result: tenfinger-decoding-accuracy reproduced. Under every reading of O1 and O2 the pooled accuracies are 71.4 to 73.0% for N-PPC (1000 trials), 68.0 to 70.0% for J-PPC and 78.0 to 78.5% for J-MC (200 trials each), within the tolerances of 70%, 66% and 75%, against a chance of 10%. The single seeds of O2 (b) range from 69.9 to 73.6%, 66.0 to 72.5% and 76.5 to 80.0%, all within tolerance.
Result: tenfinger-error-structure reproduced. Of the errors, 73.9 to 76.6% (N-PPC), 66.2 to 72.1% (J-PPC) and 62.8 to 65.3% (J-MC) confuse an adjacent finger of the same hand or the matching finger of the other hand; every value is within tolerance of 76%, 66% and 68% under every reading.
fig8b-matching-closer¶
Matching right and left finger pairs are closer in neural distance than non-matching pairs, by 1.56 units in N-PPC (95% CI 1.33 to 1.78) and 4.30 in J-MC (2.74 to 5.46), but not conclusively in J-PPC (0.27, -0.17 to 0.64).
Cross-validated squared Mahalanobis distances (rsatoolbox), per session and array, on Go-window rates; readings O1 x O3 x O4. For O4 (a) the per-partition precision is multiplied by the number of trials per condition in the partition (5), because rsatoolbox's calc_rdm_crossnobis documents that a list of noises must be "the precision of the mean from the corresponding crossvalidation fold". This is an assumption about how the paper called the toolbox; O4 (a1), the unscaled single-trial precision, is shown as a variant. The 95% CI is DABEST's (O5 a, BCa bootstrap of the difference of means, 5000 resamples, each group resampled on its own) and the percentile bootstrap (O5 b).
rdms, rdm_rows = {}, []
for o1 in ["a", "b"]:
for o3 in ["shrinkage_diag", "shrinkage_eye", "diag"]:
for o4 in ["a", "b", "a1"]:
for arr in F.ARRAYS:
M, N, R, dropped = [], [], [], []
for s in sess_of(arr):
X, y = F.array_data(s, arr, o1)
rdm, nd = F.crossnobis_rdm(X, y, o3, o4)
m, n = F.cross_hand_elements(rdm)
M += list(m); N += list(n); R.append(rdm); dropped.append(nd)
M, N = np.array(M), np.array(N)
rdms[(o1, o3, o4, arr)] = R
stat = lambda a, b, axis=-1: a.mean(axis) - b.mean(axis)
row = dict(O1=o1, O3=o3, O4=o4, array=arr, n_match=len(M), n_non=len(N), dropped_units=sum(dropped),
match_mean=M.mean(), non_mean=N.mean(), diff=N.mean() - M.mean())
for meth, lab in [("BCa", "bca"), ("percentile", "pct")]:
ci = bootstrap((N, M), stat, n_resamples=5000, method=meth, vectorized=True,
random_state=np.random.default_rng(8)).confidence_interval
row[f"{lab}_low"], row[f"{lab}_high"] = ci.low, ci.high
rdm_rows.append(row)
rdm_df = pd.DataFrame(rdm_rows)
exp = {"N-PPC": (1.56, 0.23), "J-MC": (4.30, 1.56), "J-PPC": (0.27, 0.44)}
rdm_df["in_tol"] = [abs(r["diff"] - exp[r["array"]][0]) <= exp[r["array"]][1] for _, r in rdm_df.iterrows()]
print(rdm_df.round(2).to_string())
O1 O3 O4 array n_match n_non dropped_units match_mean non_mean diff bca_low bca_high pct_low pct_high in_tol 0 a shrinkage_diag a N-PPC 50 200 1 1.07 2.49 1.42 1.19 1.64 1.19 1.64 True 1 a shrinkage_diag a J-PPC 10 40 0 1.95 2.17 0.22 -0.10 0.67 -0.13 0.62 True 2 a shrinkage_diag a J-MC 10 40 0 2.00 5.05 3.06 1.73 4.13 1.85 4.20 True 3 a shrinkage_diag b N-PPC 50 200 1 0.21 0.51 0.29 0.25 0.34 0.25 0.34 False 4 a shrinkage_diag b J-PPC 10 40 0 0.40 0.44 0.04 -0.01 0.12 -0.01 0.11 True 5 a shrinkage_diag b J-MC 10 40 0 0.43 0.94 0.51 0.27 0.70 0.28 0.71 False 6 a shrinkage_diag a1 N-PPC 50 200 1 0.21 0.50 0.28 0.24 0.33 0.24 0.33 False 7 a shrinkage_diag a1 J-PPC 10 40 0 0.39 0.43 0.04 -0.02 0.13 -0.03 0.12 True 8 a shrinkage_diag a1 J-MC 10 40 0 0.40 1.01 0.61 0.35 0.83 0.37 0.84 False 9 a shrinkage_eye a N-PPC 50 200 1 1.32 3.02 1.70 1.41 1.99 1.41 1.99 True 10 a shrinkage_eye a J-PPC 10 40 0 2.51 2.70 0.19 -0.30 0.69 -0.31 0.68 True 11 a shrinkage_eye a J-MC 10 40 0 2.19 5.19 2.99 1.32 4.27 1.46 4.39 True 12 a shrinkage_eye b N-PPC 50 200 1 0.35 0.71 0.36 0.29 0.41 0.29 0.41 False 13 a shrinkage_eye b J-PPC 10 40 0 0.65 0.71 0.05 -0.04 0.16 -0.04 0.15 True 14 a shrinkage_eye b J-MC 10 40 0 0.61 1.24 0.63 0.26 0.90 0.30 0.93 False 15 a shrinkage_eye a1 N-PPC 50 200 1 0.26 0.60 0.34 0.28 0.40 0.28 0.40 False 16 a shrinkage_eye a1 J-PPC 10 40 0 0.50 0.54 0.04 -0.06 0.14 -0.06 0.14 True 17 a shrinkage_eye a1 J-MC 10 40 0 0.44 1.04 0.60 0.26 0.85 0.29 0.88 False 18 a diag a N-PPC 50 200 1 1.07 2.49 1.42 1.19 1.64 1.19 1.64 True 19 a diag a J-PPC 10 40 0 1.95 2.17 0.22 -0.10 0.67 -0.13 0.62 True 20 a diag a J-MC 10 40 0 2.10 5.42 3.32 1.94 4.44 2.05 4.53 True 21 a diag b N-PPC 50 200 1 0.21 0.51 0.29 0.25 0.34 0.25 0.34 False 22 a diag b J-PPC 10 40 0 0.38 0.43 0.05 -0.01 0.13 -0.01 0.12 True 23 a diag b J-MC 10 40 0 0.43 1.08 0.65 0.37 0.88 0.39 0.89 False 24 a diag a1 N-PPC 50 200 1 0.21 0.50 0.28 0.24 0.33 0.24 0.33 False 25 a diag a1 J-PPC 10 40 0 0.39 0.43 0.04 -0.02 0.13 -0.03 0.12 True 26 a diag a1 J-MC 10 40 0 0.42 1.08 0.66 0.39 0.89 0.41 0.91 False 27 b shrinkage_diag a N-PPC 50 200 1 1.08 2.52 1.44 1.22 1.66 1.22 1.66 True 28 b shrinkage_diag a J-PPC 10 40 0 1.96 2.18 0.22 -0.11 0.66 -0.14 0.62 True 29 b shrinkage_diag a J-MC 10 40 0 2.00 5.13 3.13 1.74 4.21 1.88 4.30 True 30 b shrinkage_diag b N-PPC 50 200 1 0.22 0.51 0.30 0.25 0.34 0.25 0.34 False 31 b shrinkage_diag b J-PPC 10 40 0 0.40 0.44 0.04 -0.01 0.12 -0.02 0.11 True 32 b shrinkage_diag b J-MC 10 40 0 0.43 0.95 0.52 0.27 0.71 0.29 0.72 False 33 b shrinkage_diag a1 N-PPC 50 200 1 0.22 0.50 0.29 0.24 0.33 0.24 0.33 False 34 b shrinkage_diag a1 J-PPC 10 40 0 0.39 0.44 0.04 -0.02 0.13 -0.03 0.12 True 35 b shrinkage_diag a1 J-MC 10 40 0 0.40 1.03 0.63 0.35 0.84 0.38 0.86 False 36 b shrinkage_eye a N-PPC 50 200 1 1.33 3.03 1.71 1.42 1.99 1.41 1.99 True 37 b shrinkage_eye a J-PPC 10 40 0 2.52 2.70 0.19 -0.31 0.69 -0.31 0.69 True 38 b shrinkage_eye a J-MC 10 40 0 2.20 5.21 3.01 1.33 4.30 1.46 4.42 True 39 b shrinkage_eye b N-PPC 50 200 1 0.35 0.71 0.36 0.29 0.42 0.29 0.42 False 40 b shrinkage_eye b J-PPC 10 40 0 0.65 0.71 0.05 -0.04 0.16 -0.04 0.15 True 41 b shrinkage_eye b J-MC 10 40 0 0.61 1.24 0.63 0.26 0.91 0.30 0.94 False 42 b shrinkage_eye a1 N-PPC 50 200 1 0.27 0.61 0.34 0.28 0.40 0.28 0.40 False 43 b shrinkage_eye a1 J-PPC 10 40 0 0.50 0.54 0.04 -0.06 0.14 -0.06 0.14 True 44 b shrinkage_eye a1 J-MC 10 40 0 0.44 1.04 0.60 0.27 0.86 0.29 0.88 False 45 b diag a N-PPC 50 200 1 1.08 2.52 1.44 1.22 1.66 1.22 1.66 True 46 b diag a J-PPC 10 40 0 1.96 2.18 0.22 -0.11 0.66 -0.14 0.62 True 47 b diag a J-MC 10 40 0 2.11 5.49 3.38 1.95 4.53 2.07 4.61 True 48 b diag b N-PPC 50 200 1 0.22 0.51 0.30 0.25 0.34 0.25 0.34 False 49 b diag b J-PPC 10 40 0 0.38 0.44 0.05 -0.01 0.13 -0.01 0.12 True 50 b diag b J-MC 10 40 0 0.43 1.09 0.66 0.37 0.89 0.39 0.91 False 51 b diag a1 N-PPC 50 200 1 0.22 0.50 0.29 0.24 0.33 0.24 0.33 False 52 b diag a1 J-PPC 10 40 0 0.39 0.44 0.04 -0.02 0.13 -0.03 0.12 True 53 b diag a1 J-MC 10 40 0 0.42 1.10 0.68 0.39 0.91 0.41 0.92 False
def verdict(g):
d = g.set_index("array")
pos = all(d.loc[a, "diff"] > 0 and d.loc[a, "bca_low"] > 0 for a in ["N-PPC", "J-MC"])
jppc = d.loc["J-PPC", "bca_low"] <= 0 <= d.loc["J-PPC", "bca_high"]
tol = bool(d["in_tol"].all())
return pd.Series(dict(NPPC_JMC_positive_CI_excl_0=pos, JPPC_CI_includes_0=jppc, all_in_tol=tol,
N_PPC=round(d.loc["N-PPC", "diff"], 2), J_MC=round(d.loc["J-MC", "diff"], 2), J_PPC=round(d.loc["J-PPC", "diff"], 2)))
print("per reading (CI by O5 a, BCa):")
print(rdm_df.groupby(["O1", "O3", "O4"]).apply(verdict).to_string())
jp = rdm_df[rdm_df.array == "J-PPC"]
print("J-PPC percentile CI includes 0 in", int(((jp.pct_low <= 0) & (jp.pct_high >= 0)).sum()), "of", len(jp), "readings")
key = ("a", "shrinkage_diag", "a")
fig, axes = plt.subplots(1, 2, figsize=(11, 4))
A = np.mean(rdms[key + ("N-PPC",)], axis=0)
im = axes[0].imshow(A, cmap="Blues")
axes[0].set_xticks(range(10), F.CONDITIONS, fontsize=8); axes[0].set_yticks(range(10), F.CONDITIONS, fontsize=8)
axes[0].set_title("N-PPC distances, mean of 10 sessions (O1 a, O3 a, O4 a)")
fig.colorbar(im, ax=axes[0], shrink=0.8)
off = ~np.eye(10, dtype=bool)
print(f"N-PPC mean RDM range, off-diagonal elements {A[off].min():.2f} to {A[off].max():.2f}")
for i, arr in enumerate(F.ARRAYS):
M, N = [], []
for r in rdms[key + (arr,)]:
m, n = F.cross_hand_elements(r); M += list(m); N += list(n)
jit = np.random.default_rng(9)
axes[1].scatter(i - 0.15 + 0.05 * jit.standard_normal(len(M)), M, s=12, color=BLUE, label="matching" if i == 0 else None)
axes[1].scatter(i + 0.15 + 0.05 * jit.standard_normal(len(N)), N, s=12, color=ORANGE, label="non-matching" if i == 0 else None)
print(f"{arr}: matching {min(M):.2f} to {max(M):.2f}, non-matching {min(N):.2f} to {max(N):.2f}")
axes[1].set_xticks(range(3), list(F.ARRAYS)); axes[1].set_ylabel("cross-validated distance (per neuron)")
axes[1].legend(frameon=False); axes[1].set_title("different-hand RDM elements, all sessions")
plt.show()
per reading (CI by O5 a, BCa):
NPPC_JMC_positive_CI_excl_0 JPPC_CI_includes_0 all_in_tol N_PPC J_MC J_PPC
O1 O3 O4
a diag a True True True 1.42 3.32 0.22
a1 True True False 0.28 0.66 0.04
b True True False 0.29 0.65 0.05
shrinkage_diag a True True True 1.42 3.06 0.22
a1 True True False 0.28 0.61 0.04
b True True False 0.29 0.51 0.04
shrinkage_eye a True True True 1.70 2.99 0.19
a1 True True False 0.34 0.60 0.04
b True True False 0.36 0.63 0.05
b diag a True True True 1.44 3.38 0.22
a1 True True False 0.29 0.68 0.04
b True True False 0.30 0.66 0.05
shrinkage_diag a True True True 1.44 3.13 0.22
a1 True True False 0.29 0.63 0.04
b True True False 0.30 0.52 0.04
shrinkage_eye a True True True 1.71 3.01 0.19
a1 True True False 0.34 0.60 0.04
b True True False 0.36 0.63 0.05
J-PPC percentile CI includes 0 in 18 of 18 readings
N-PPC mean RDM range, off-diagonal elements 0.57 to 3.90
N-PPC: matching 0.06 to 2.36, non-matching -0.07 to 6.19
J-PPC: matching 0.76 to 2.60, non-matching 1.11 to 3.75
J-MC: matching 0.63 to 5.48, non-matching 1.29 to 11.47
# The scale of the distances under O4 (a), (a1) and (b): N-PPC RDM range, mean over sessions and per session (O1 a, O3 a)
rng_rows = []
for o4 in ["a", "a1", "b"]:
R = rdms[("a", "shrinkage_diag", o4, "N-PPC")]
off = ~np.eye(10, dtype=bool)
A = np.mean(R, axis=0)
rng_rows.append(dict(O4=o4, mean_rdm_min=A[off].min(), mean_rdm_max=A[off].max(),
per_session_min=[round(r[off].min(), 2) for r in R], per_session_max=[round(r[off].max(), 2) for r in R]))
print(pd.DataFrame(rng_rows).round(2).to_string())
# Synthetic check of the scale: 10 conditions x 10 trials, 100 units with unit-variance Gaussian noise and known
# condition means; ratio of the mean estimated distance to the true squared Mahalanobis distance per neuron
# (in single-trial noise units), noise estimator 'diag'.
syn = np.random.default_rng(0)
mu = syn.normal(0, 0.5, (10, 100))
y_syn = np.repeat(F.CONDITIONS, 10)
X_syn = mu[np.repeat(np.arange(10), 10)] + syn.normal(0, 1, (100, 100))
true = ((mu[:, None] - mu[None]) ** 2).sum(-1) / 100
iu = np.triu_indices(10, 1)
print("synthetic: estimated / true distance:",
{o4: round(float(np.mean(F.crossnobis_rdm(X_syn, y_syn, "diag", o4)[0][iu]) / np.mean(true[iu])), 2) for o4 in ["a", "a1", "b"]})
O4 mean_rdm_min mean_rdm_max per_session_min per_session_max
0 a 0.57 3.90 [0.21, -0.07, 0.46, 0.06, 0.29, 0.37, 0.16, 0.49, 0.62, 0.17] [5.89, 5.57, 6.19, 5.63, 2.71, 3.11, 3.56, 4.45, 4.38, 3.29]
1 a1 0.11 0.78 [0.04, -0.01, 0.09, 0.01, 0.06, 0.07, 0.03, 0.1, 0.12, 0.03] [1.18, 1.11, 1.24, 1.13, 0.54, 0.62, 0.71, 0.89, 0.88, 0.66]
2 b 0.12 0.79 [0.07, 0.03, 0.05, 0.05, 0.1, 0.09, 0.05, 0.12, 0.11, 0.06] [1.23, 1.14, 1.26, 1.11, 0.57, 0.57, 0.72, 0.87, 0.86, 0.75]
synthetic: estimated / true distance: {'a': 5.17, 'a1': 1.03, 'b': 1.03}
Result: reproduced, resting on an assumption. Under O4 (a), with the toolbox's per-partition precision of the mean, the non-matching minus matching difference is 1.42 to 1.71 in N-PPC, 2.99 to 3.38 in J-MC and 0.19 to 0.22 in J-PPC over the readings of O1 and O3; N-PPC and J-MC have positive differences with BCa CIs above 0, the J-PPC CI includes 0, and all values are within tolerance (the largest N-PPC value, 1.71 under shrinkage toward identity, against an upper bound of 1.79). Under O4 (b) and the variant O4 (a1) the same pattern holds (positive CIs excluding 0 for N-PPC and J-MC, J-PPC CI including 0) but the differences are about a fifth as large (N-PPC 0.28 to 0.36, J-MC 0.51 to 0.68), outside tolerance, which would be partially reproduced. The outcome therefore rests on how rsatoolbox's noise argument was scaled, which the paper does not state. The cell above shows the scale directly: on synthetic data with a known distance, O4 (a) returns 5.17 times the distance in single-trial noise units, O4 (a1) and (b) 1.03 times; the N-PPC mean RDM spans 0.57 to 3.90 under O4 (a) and 0.11 to 0.78 under O4 (a1) and 0.12 to 0.79 under O4 (b). The text of the paper gives no colour-bar range for Figure 8a, so which scale its figure uses cannot be checked here.
fig8a-npcc-index-distinct¶
In N-PPC the inter-finger geometry is similar across hands, with the index finger relatively distinct from the middle and ring fingers in both hands.
From the N-PPC distances of fig8b-matching-closer averaged over the 10 sessions; reproduced when in each hand the index-middle and index-ring distances both exceed the middle-ring distance, under one reading. MDS with Generalized Procrustes alignment, as the Methods describe, is drawn for index, middle and ring.
ix = {c: i for i, c in enumerate(F.CONDITIONS)}
rows = []
for (o1, o3, o4, arr), R in rdms.items():
if arr != "N-PPC":
continue
A = np.mean(R, axis=0)
r = dict(O1=o1, O3=o3, O4=o4)
ok = True
for h in "RL":
im_, ir_, mr_ = A[ix[h + "i"], ix[h + "m"]], A[ix[h + "i"], ix[h + "r"]], A[ix[h + "m"], ix[h + "r"]]
r.update({f"{h}_IM": im_, f"{h}_IR": ir_, f"{h}_MR": mr_})
ok &= (im_ > mr_) and (ir_ > mr_)
# per-session count of sessions where the pattern holds in both hands
n_ok = 0
for Rs in R:
n_ok += all(Rs[ix[h + "i"], ix[h + "m"]] > Rs[ix[h + "m"], ix[h + "r"]] and Rs[ix[h + "i"], ix[h + "r"]] > Rs[ix[h + "m"], ix[h + "r"]] for h in "RL")
r.update(pattern_both_hands=ok, sessions_with_pattern=n_ok)
rows.append(r)
print(pd.DataFrame(rows).round(2).to_string())
sub = [c for c in F.CONDITIONS if c[1] in "imr"]
si = [ix[c] for c in sub]
confs = [F.classical_mds(Rs[np.ix_(si, si)]) for Rs in rdms[("a", "shrinkage_diag", "a", "N-PPC")]]
aligned = np.array(F.gpa(confs))
mean, se = aligned.mean(axis=0), aligned.std(axis=0, ddof=1) / np.sqrt(len(aligned))
fig, ax = plt.subplots(figsize=(4.5, 4.5))
for k, c in enumerate(sub):
col = BLUE if c[0] == "R" else ORANGE
ax.errorbar(mean[k, 0], mean[k, 1], xerr=se[k, 0], yerr=se[k, 1], fmt="o", color=col, ms=8)
ax.annotate(c, mean[k], textcoords="offset points", xytext=(6, 6))
for h, col in [("R", BLUE), ("L", ORANGE)]:
pts = np.array([mean[sub.index(h + f)] for f in "imr"])
ax.plot(pts[:, 0], pts[:, 1], color=col, lw=2)
ax.set_title("N-PPC MDS, GPA over 10 sessions (mean, SE)")
ax.set_xticks([]); ax.set_yticks([])
plt.show()
print("MDS coordinates (mean over sessions, arbitrary scale):")
print(pd.DataFrame(mean, index=sub, columns=["x", "y"]).round(3).to_string())
O1 O3 O4 R_IM R_IR R_MR L_IM L_IR L_MR pattern_both_hands sessions_with_pattern 0 a shrinkage_diag a 2.16 2.96 0.62 1.30 2.23 0.57 True 8 1 a shrinkage_diag b 0.43 0.58 0.13 0.26 0.45 0.12 True 10 2 a shrinkage_diag a1 0.43 0.59 0.12 0.26 0.45 0.11 True 8 3 a shrinkage_eye a 2.52 3.50 0.75 1.54 2.67 0.79 True 7 4 a shrinkage_eye b 0.58 0.78 0.24 0.40 0.65 0.23 True 8 5 a shrinkage_eye a1 0.50 0.70 0.15 0.31 0.53 0.16 True 7 6 a diag a 2.16 2.96 0.62 1.30 2.23 0.57 True 8 7 a diag b 0.43 0.58 0.13 0.26 0.45 0.12 True 10 8 a diag a1 0.43 0.59 0.12 0.26 0.45 0.11 True 8 9 b shrinkage_diag a 2.18 2.98 0.63 1.34 2.27 0.59 True 9 10 b shrinkage_diag b 0.44 0.59 0.13 0.27 0.47 0.12 True 10 11 b shrinkage_diag a1 0.44 0.60 0.13 0.27 0.45 0.12 True 9 12 b shrinkage_eye a 2.52 3.51 0.76 1.54 2.68 0.80 True 7 13 b shrinkage_eye b 0.58 0.78 0.24 0.40 0.65 0.23 True 8 14 b shrinkage_eye a1 0.50 0.70 0.15 0.31 0.54 0.16 True 7 15 b diag a 2.18 2.98 0.63 1.34 2.27 0.59 True 9 16 b diag b 0.44 0.59 0.13 0.27 0.47 0.12 True 10 17 b diag a1 0.44 0.60 0.13 0.27 0.45 0.12 True 9
MDS coordinates (mean over sessions, arbitrary scale):
x y
Ri -0.468 0.130
Rm 0.137 0.206
Rr 0.295 0.253
Li -0.341 -0.218
Lm 0.084 -0.193
Lr 0.294 -0.178
Result: reproduced. In the N-PPC distances averaged over the 10 sessions, under every reading of O1, O3 and O4, the index-middle and index-ring distances exceed the middle-ring distance in both hands (for example, O1 (a), O4 (a) with diagonal shrinkage: right hand 2.16 and 2.96 against 0.62, left hand 1.30 and 2.23 against 0.57). In single sessions the pattern holds in both hands in 7 to 10 of the 10 sessions, depending on the reading.
fig8e-hand-ccgp, fig8e-finger-ccgp and jppc-finger-ccgp-low¶
Hand (left against right) decoders trained on four finger types generalize to the held-out finger type significantly above chance (p < 0.001, permutation test), with hand CCGP close to standard decoding accuracy, in N-PPC, J-MC and J-PPC.
Finger-type decoders trained on one hand generalize to the other hand significantly above chance (p < 0.001), with finger-type CCGP close to standard decoding accuracy, in N-PPC and J-MC.
J-PPC finger representations are less compositional; finger-type CCGP (36%) is above chance (p < 0.001) but much lower than within-hand accuracy (65%).
Pseudo-population of the array's neurons over sessions (O6), Go-window rates, decoders O7 (a, diagonal LDA; b, linear SVM), 1001 label shuffles, p = (1 + shuffles reaching the observed CCGP) / 1002, so p < 0.001 means no shuffle reached it. Standard accuracy (8-fold, O2 a) is computed on two populations, each under every O1 x O6 x O7: (1) the same pseudo-population as CCGP, which this attempt adopts, since the Methods (Compositionality of finger representations) call it "an upper bound on CCGP", a bound that holds only for the same neurons, and the CCGP literature it cites (ref 60) computes both on one pseudo-population; (2) per session (Methods, Offline classification: "a separate linear classifier for each session", trials pooled), which does not depend on O6. The paper does not say which it used, so the choice is an assumption. For O6 (b) the observed CCGP is the mean over 100 matchings while each null value uses one matching, which widens the null relative to the observed statistic and makes p conservative; it is stated as an assumption. The number of null values reaching the observed CCGP is printed. Closeness by O9: (a) CCGP at least 0.9 of standard accuracy; (b) ratio above 0.55. Finger standard accuracy by O8: (a) within hand, (b) over both hands.
t = time.time()
ccgp_rows = []
for o1 in ["a", "b"]:
for o6 in ["a", "b"]:
for o7 in ["lda", "svm"]:
for k, arr in enumerate(F.ARRAYS):
Xs, ys = zip(*[F.array_data(s, arr, o1) for s in sess_of(arr)])
rng = np.random.default_rng(1000 + k) # this section's own generator, per array
r = F.run_ccgp(list(Xs), list(ys), o6, o7, rng)
ccgp_rows.append(dict(O1=o1, O6=o6, O7=o7, array=arr, **r))
ccgp_df = pd.DataFrame(ccgp_rows)
ccgp_df["hand_ratio"] = ccgp_df.hand / ccgp_df.std_hand
ccgp_df["finger_ratio_O8a"] = ccgp_df.finger / ccgp_df.std_within
ccgp_df["finger_ratio_O8b"] = ccgp_df.finger / ccgp_df.std_both
print(f"{time.time() - t:.0f} s")
print(ccgp_df.round(3).to_string())
749 s O1 O6 O7 array hand finger std_hand std_within std_both p_hand p_finger n_null_ge_hand n_null_ge_finger null_hand_mean null_hand_max null_finger_mean null_finger_max n_neurons hand_ratio finger_ratio_O8a finger_ratio_O8b 0 a a lda N-PPC 0.980 0.940 1.000 0.980 0.980 0.001 0.001 0 0 0.500 0.68 0.201 0.37 943 0.980 0.959 0.959 1 a a lda J-PPC 0.980 0.290 0.980 0.870 0.790 0.001 0.038 0 37 0.499 0.67 0.201 0.36 315 1.000 0.333 0.367 2 a a lda J-MC 0.870 0.780 0.890 0.980 0.910 0.001 0.001 0 0 0.498 0.67 0.200 0.35 252 0.978 0.796 0.857 3 a a svm N-PPC 0.980 0.930 1.000 0.960 0.970 0.001 0.001 0 0 0.501 0.69 0.200 0.35 943 0.980 0.969 0.959 4 a a svm J-PPC 0.980 0.360 0.980 0.840 0.760 0.001 0.002 0 1 0.500 0.69 0.201 0.38 315 1.000 0.429 0.474 5 a a svm J-MC 0.860 0.790 0.930 0.960 0.930 0.001 0.001 0 0 0.500 0.70 0.199 0.38 252 0.925 0.823 0.849 6 a b lda N-PPC 0.987 0.947 0.997 0.994 0.992 0.001 0.001 0 0 0.498 0.68 0.196 0.33 943 0.990 0.953 0.955 7 a b lda J-PPC 0.984 0.332 0.986 0.845 0.773 0.001 0.003 0 2 0.502 0.67 0.200 0.40 315 0.998 0.393 0.430 8 a b lda J-MC 0.862 0.781 0.886 0.955 0.922 0.001 0.001 0 0 0.501 0.68 0.199 0.38 252 0.973 0.818 0.847 9 a b svm N-PPC 0.974 0.918 0.997 0.993 0.990 0.001 0.001 0 0 0.500 0.72 0.197 0.33 943 0.976 0.925 0.927 10 a b svm J-PPC 0.970 0.337 0.977 0.828 0.751 0.001 0.004 0 3 0.499 0.69 0.199 0.37 315 0.993 0.407 0.449 11 a b svm J-MC 0.865 0.735 0.943 0.954 0.949 0.001 0.001 0 0 0.500 0.68 0.199 0.37 252 0.917 0.770 0.775 12 b a lda N-PPC 0.990 0.950 1.000 0.980 0.980 0.001 0.001 0 0 0.500 0.67 0.201 0.38 930 0.990 0.969 0.969 13 b a lda J-PPC 0.990 0.320 0.980 0.890 0.800 0.001 0.016 0 15 0.499 0.68 0.201 0.38 312 1.010 0.360 0.400 14 b a lda J-MC 0.860 0.780 0.900 0.970 0.910 0.001 0.001 0 0 0.498 0.67 0.200 0.35 250 0.956 0.804 0.857 15 b a svm N-PPC 0.980 0.940 1.000 0.970 0.970 0.001 0.001 0 0 0.502 0.69 0.199 0.35 930 0.980 0.969 0.969 16 b a svm J-PPC 0.980 0.380 0.990 0.850 0.760 0.001 0.001 0 0 0.500 0.70 0.201 0.37 312 0.990 0.447 0.500 17 b a svm J-MC 0.860 0.790 0.930 0.960 0.930 0.001 0.001 0 0 0.498 0.71 0.198 0.37 250 0.925 0.823 0.849 18 b b lda N-PPC 0.987 0.954 0.997 0.995 0.994 0.001 0.001 0 0 0.497 0.68 0.196 0.32 930 0.989 0.960 0.961 19 b b lda J-PPC 0.982 0.339 0.985 0.845 0.771 0.001 0.002 0 1 0.501 0.69 0.201 0.39 312 0.997 0.401 0.440 20 b b lda J-MC 0.861 0.783 0.886 0.955 0.924 0.001 0.001 0 0 0.500 0.68 0.199 0.38 250 0.972 0.820 0.848 21 b b svm N-PPC 0.976 0.927 0.998 0.992 0.992 0.001 0.001 0 0 0.498 0.72 0.197 0.33 930 0.978 0.935 0.935 22 b b svm J-PPC 0.971 0.358 0.978 0.828 0.749 0.001 0.002 0 1 0.500 0.69 0.200 0.38 312 0.993 0.432 0.478 23 b b svm J-MC 0.860 0.745 0.943 0.953 0.948 0.001 0.001 0 0 0.499 0.70 0.199 0.37 250 0.913 0.781 0.786
# Standard accuracy computed per session (a classifier per session, stratified 8-fold, O2 a, trials pooled
# over sessions; Methods, Offline classification), under every O1 x O7. It does not depend on O6, so the same
# values join every O6 row of ccgp_df.
ps_rows = []
for o1 in ["a", "b"]:
for o7 in ["lda", "svm"]:
for arr in F.ARRAYS:
r = {}
for tgt in ["hand", "finger_within", "finger_both"]:
a, n = [], []
for s in sess_of(arr):
X, y = F.array_data(s, arr, o1)
a.append(F.standard_accuracy(X, y, tgt, o7)); n.append(len(y))
r["ps_" + tgt] = np.average(a, weights=n)
ps_rows.append(dict(O1=o1, O7=o7, array=arr, **r))
ccgp_df = ccgp_df.drop(columns=[c for c in ccgp_df if c.startswith("ps_")]).merge(pd.DataFrame(ps_rows), on=["O1", "O7", "array"])
# Each population of the standard accuracy: 'pseudo' (the same pseudo-population as CCGP, assumed) and 'per_session'
STD = {"pseudo": ("std_hand", "std_within", "std_both"), "per_session": ("ps_hand", "ps_finger_within", "ps_finger_both")}
print("number of the 1001 null values reaching the observed CCGP (p = (1 + n) / 1002):")
print(ccgp_df[["O1", "O6", "O7", "array", "n_null_ge_hand", "p_hand", "n_null_ge_finger", "p_finger"]].to_string())
for pop, (sh, sw, sb) in STD.items():
d = ccgp_df.assign(hr=ccgp_df.hand / ccgp_df[sh], fa=ccgp_df.finger / ccgp_df[sw], fb=ccgp_df.finger / ccgp_df[sb])
g = d.groupby(["O1", "O6", "O7"])
print(f"\n[standard accuracy: {pop}] fig8e-hand-ccgp: every array p < 0.001, hand ratio >= 0.9 (O9 a) / > 0.55 (O9 b):")
print(g.apply(lambda x: pd.Series(dict(all_p=bool((x.p_hand < 0.001).all()), O9a=bool((x.hr >= 0.9).all()),
O9b=bool((x.hr > 0.55).all()), min_ratio=round(x.hr.min(), 3)))).to_string())
nj = d[d.array.isin(["N-PPC", "J-MC"])].groupby(["O1", "O6", "O7"])
print(f"\n[standard accuracy: {pop}] fig8e-finger-ccgp (N-PPC and J-MC): p < 0.001 and finger ratio by O8 a/b and O9 a/b:")
print(nj.apply(lambda x: pd.Series(dict(all_p=bool((x.p_finger < 0.001).all()),
O8a_O9a=bool((x.fa >= 0.9).all()), O8b_O9a=bool((x.fb >= 0.9).all()),
O8a_O9b=bool((x.fa > 0.55).all()), O8b_O9b=bool((x.fb > 0.55).all()),
NPPC_ratio_O8a=round(x.set_index("array").loc["N-PPC", "fa"], 3),
JMC_ratio_O8a=round(x.set_index("array").loc["J-MC", "fa"], 3)))).to_string())
jp = ccgp_df[ccgp_df.array == "J-PPC"][["O1", "O6", "O7", "finger", "p_finger", "n_null_ge_finger",
"std_within", "std_both", "ps_finger_within", "ps_finger_both"]].copy()
jp["ccgp_in_tol"] = (jp.finger * 100 - 36).abs() <= 3.6
for pop, col in [("pseudo", "std_within"), ("per_session", "ps_finger_within")]:
within_ok = (jp[col] * 100 - 65).abs() <= 6.5
low = jp.finger < 0.9 * jp[col]
jp[f"{pop}_within_in_tol"] = within_ok
jp[f"{pop}_verdict"] = np.where((jp.p_finger < 0.001) & jp.ccgp_in_tol & within_ok, "reproduced",
np.where((jp.p_finger < 0.001) & low, "partially", "no"))
print("\njppc-finger-ccgp-low (expected CCGP 36 +/- 3.6%, within-hand 65 +/- 6.5%), per O1 x O6 x O7 and population of the within-hand accuracy:")
print(jp.round(3).to_string())
number of the 1001 null values reaching the observed CCGP (p = (1 + n) / 1002):
O1 O6 O7 array n_null_ge_hand p_hand n_null_ge_finger p_finger
0 a a lda N-PPC 0 0.000998 0 0.000998
1 a a lda J-PPC 0 0.000998 37 0.037924
2 a a lda J-MC 0 0.000998 0 0.000998
3 a a svm N-PPC 0 0.000998 0 0.000998
4 a a svm J-PPC 0 0.000998 1 0.001996
5 a a svm J-MC 0 0.000998 0 0.000998
6 a b lda N-PPC 0 0.000998 0 0.000998
7 a b lda J-PPC 0 0.000998 2 0.002994
8 a b lda J-MC 0 0.000998 0 0.000998
9 a b svm N-PPC 0 0.000998 0 0.000998
10 a b svm J-PPC 0 0.000998 3 0.003992
11 a b svm J-MC 0 0.000998 0 0.000998
12 b a lda N-PPC 0 0.000998 0 0.000998
13 b a lda J-PPC 0 0.000998 15 0.015968
14 b a lda J-MC 0 0.000998 0 0.000998
15 b a svm N-PPC 0 0.000998 0 0.000998
16 b a svm J-PPC 0 0.000998 0 0.000998
17 b a svm J-MC 0 0.000998 0 0.000998
18 b b lda N-PPC 0 0.000998 0 0.000998
19 b b lda J-PPC 0 0.000998 1 0.001996
20 b b lda J-MC 0 0.000998 0 0.000998
21 b b svm N-PPC 0 0.000998 0 0.000998
22 b b svm J-PPC 0 0.000998 1 0.001996
23 b b svm J-MC 0 0.000998 0 0.000998
[standard accuracy: pseudo] fig8e-hand-ccgp: every array p < 0.001, hand ratio >= 0.9 (O9 a) / > 0.55 (O9 b):
all_p O9a O9b min_ratio
O1 O6 O7
a a lda True True True 0.978
svm True True True 0.925
b lda True True True 0.973
svm True True True 0.917
b a lda True True True 0.956
svm True True True 0.925
b lda True True True 0.972
svm True True True 0.913
[standard accuracy: pseudo] fig8e-finger-ccgp (N-PPC and J-MC): p < 0.001 and finger ratio by O8 a/b and O9 a/b:
all_p O8a_O9a O8b_O9a O8a_O9b O8b_O9b NPPC_ratio_O8a JMC_ratio_O8a
O1 O6 O7
a a lda True False False True True 0.959 0.796
svm True False False True True 0.969 0.823
b lda True False False True True 0.953 0.818
svm True False False True True 0.925 0.770
b a lda True False False True True 0.969 0.804
svm True False False True True 0.969 0.823
b lda True False False True True 0.960 0.820
svm True False False True True 0.935 0.781
[standard accuracy: per_session] fig8e-hand-ccgp: every array p < 0.001, hand ratio >= 0.9 (O9 a) / > 0.55 (O9 b):
all_p O9a O9b min_ratio
O1 O6 O7
a a lda True True True 1.012
svm True True True 0.940
b lda True True True 1.003
svm True True True 0.945
b a lda True True True 0.994
svm True True True 0.945
b lda True True True 0.995
svm True True True 0.946
[standard accuracy: per_session] fig8e-finger-ccgp (N-PPC and J-MC): p < 0.001 and finger ratio by O8 a/b and O9 a/b:
all_p O8a_O9a O8b_O9a O8a_O9b O8b_O9b NPPC_ratio_O8a JMC_ratio_O8a
O1 O6 O7
a a lda True False True True True 1.146 0.897
svm True True True True True 1.299 0.908
b lda True False True True True 1.155 0.898
svm True False False True True 1.283 0.845
b a lda True False True True True 1.154 0.897
svm True True True True True 1.309 0.903
b lda True False True True True 1.160 0.900
svm True False False True True 1.291 0.851
jppc-finger-ccgp-low (expected CCGP 36 +/- 3.6%, within-hand 65 +/- 6.5%), per O1 x O6 x O7 and population of the within-hand accuracy:
O1 O6 O7 finger p_finger n_null_ge_finger std_within std_both ps_finger_within ps_finger_both ccgp_in_tol pseudo_within_in_tol pseudo_verdict per_session_within_in_tol per_session_verdict
1 a a lda 0.290 0.038 37 0.870 0.790 0.690 0.695 False False no True no
4 a a svm 0.360 0.002 1 0.840 0.760 0.685 0.530 True False no True no
7 a b lda 0.332 0.003 2 0.845 0.773 0.690 0.695 True False no True no
10 a b svm 0.337 0.004 3 0.828 0.751 0.685 0.530 True False no True no
13 b a lda 0.320 0.016 15 0.890 0.800 0.685 0.700 False False no True no
16 b a svm 0.380 0.001 0 0.850 0.760 0.680 0.480 True False partially True reproduced
19 b b lda 0.339 0.002 1 0.845 0.771 0.685 0.700 True False no True no
22 b b svm 0.358 0.002 1 0.828 0.749 0.680 0.480 True False no True no
d = ccgp_df[(ccgp_df.O1 == "a") & (ccgp_df.O6 == "a") & (ccgp_df.O7 == "lda")].set_index("array")
fig, axes = plt.subplots(1, 2, figsize=(10, 3.8), sharey=True)
for ax, (cc, sd, nul, title) in zip(axes, [("hand", "std_hand", "null_hand_mean", "hand"), ("finger", "std_within", "null_finger_mean", "finger type")]):
x = np.arange(3)
ax.scatter(x - 0.1, d[cc], color=BLUE, s=60, label="CCGP")
ax.scatter(x + 0.1, d[sd], color=ORANGE, marker="^", s=60, label="standard accuracy")
ax.scatter(x, d[nul], color="gray", marker="_", s=200, label="shuffle mean")
ax.set_xticks(x, d.index); ax.set_ylim(0, 1.05); ax.set_title(f"{title} (O1 a, O6 a, O7 a)")
axes[0].set_ylabel("accuracy"); axes[0].legend(frameon=False, fontsize=8)
plt.show()
print("all plotted values lie in [0, 1]")
all plotted values lie in [0, 1]
Result: fig8e-hand-ccgp reproduced. In all three arrays, under every reading of O1, O6 and O7, no shuffle reached the observed hand CCGP (0 of 1001, p = 1/1002 < 0.001), and hand CCGP is at least 0.9 of the hand accuracy on the same pseudo-population (O9 a; the smallest ratio is 0.913), so it is also above 0.55 (O9 b). With the per-session standard accuracy the smallest ratio is 0.940, so the outcome is the same under both populations. Hand CCGP is 0.970 to 0.990 in N-PPC and J-PPC and 0.860 to 0.870 in J-MC, against shuffle means of about 0.50.
Result: fig8e-finger-ccgp reproduced under O9 (b). In N-PPC and J-MC no shuffle reached the finger-type CCGP (0 of 1001, p < 0.001) under every reading. N-PPC finger CCGP (0.918 to 0.954) is at least 0.9 of its standard accuracy on the pseudo-population, but J-MC finger CCGP (0.735 to 0.790) is 0.770 to 0.823 of the within-hand accuracy (O8 a) and 0.775 to 0.857 of the accuracy over both hands (O8 b), below 0.9 in both; against the ratio of O9 (b), 0.55, both arrays pass. With the per-session standard accuracy both arrays also pass O9 (a) in some readings (for example O1 a, O6 a, SVM, J-MC ratio 0.908 by O8 a), and N-PPC CCGP exceeds the per-session accuracy (ratio 1.15 to 1.31), which a true upper bound would not allow. The criterion requires closeness under one reading, so the claim is reproduced under both populations.
Result: jppc-finger-ccgp-low partially reproduced, resting on an assumption. J-PPC finger-type CCGP is 0.290 to 0.380 against a shuffle mean of 0.20, and far below the within-hand accuracy on the pseudo-population (0.828 to 0.890; ratio 0.333 to 0.447). The permutation p is below 0.001 in only one of the eight readings (O1 b, O6 a, linear SVM: 0 of 1001 null values reach the CCGP of 0.380, p = 1/1002); in the others 1 to 37 null values reach it (p = 0.002 to 0.038), so the paper's p < 0.001 holds in one reading, at the floor of the test, and is fragile. In that reading the CCGP is within tolerance of 36% but the within-hand accuracy on the pseudo-population (0.850) is outside the tolerance of 65%, so by the claim's rule (p < 0.001 and CCGP below 0.9 of the within-hand accuracy, a value out of tolerance) it is partially reproduced. With the per-session within-hand accuracy, now run under every O1 x O7 (0.680 to 0.690), the value is within tolerance in all eight readings, and the same reading (O1 b, O6 a, SVM) is reproduced. The two populations thus give different outcomes. This attempt keeps the pseudo-population, for the reasons above (the Methods' upper bound, ref 60, and N-PPC CCGP above the per-session accuracy), so the outcome rests on that assumption.
figs11-jppc-lower-matched-size¶
Even at matched neural population size, J-PPC finger-type CCGP is lower than that of N-PPC and J-MC.
Finger-type CCGP for 100 random subsets of the pooled neurons at sizes 50 to 250 in steps of 50 (up to the smallest pooled population, J-MC, 252 neurons under O1 a), under the readings O1 x O6 x O7 of fig8e-finger-ccgp (the pooled J-MC population is 250 neurons under O1 b, so 250 remains the largest size).
t = time.time()
sizes = [50, 100, 150, 200, 250]
ndc_rows = []
for o1 in ["a", "b"]:
for o6 in ["a", "b"]:
for o7 in ["lda", "svm"]:
for k, arr in enumerate(F.ARRAYS):
Xs, ys = zip(*[F.array_data(s, arr, o1) for s in sess_of(arr)])
rng = np.random.default_rng(2000 + k)
res = F.ndc_finger(list(Xs), list(ys), sizes, 100, o6, o7, rng)
for n, (m, sd) in res.items():
ndc_rows.append(dict(O1=o1, O6=o6, O7=o7, array=arr, size=n, mean=m, sd=sd))
ndc = pd.DataFrame(ndc_rows)
print(f"{time.time() - t:.0f} s")
wide = ndc.pivot_table(index=["O1", "O6", "O7", "size"], columns="array", values="mean")
wide["JPPC_lowest"] = (wide["J-PPC"] < wide["N-PPC"]) & (wide["J-PPC"] < wide["J-MC"])
print(wide.round(3).to_string())
print("J-PPC lowest at every size, per reading:")
print(wide.groupby(level=[0, 1, 2]).JPPC_lowest.all().to_string())
fig, ax = plt.subplots(figsize=(5.5, 3.8))
for arr, col, ls in [("N-PPC", BLUE, "-"), ("J-MC", ORANGE, "-"), ("J-PPC", "gray", "-")]:
q = ndc[(ndc.O1 == "a") & (ndc.O6 == "a") & (ndc.O7 == "lda") & (ndc.array == arr)]
ax.errorbar(q["size"], q["mean"], yerr=q["sd"], color=col, lw=2, marker="o", label=arr)
ax.axhline(0.2, color="gray", lw=1, ls=":"); ax.set_ylim(0, 1)
ax.set_xlabel("neurons"); ax.set_ylabel("finger-type CCGP"); ax.legend(frameon=False); ax.set_title("O1 a, O6 a, O7 a (mean, SD of 100 subsets)")
plt.show()
77 s
array J-MC J-PPC N-PPC JPPC_lowest
O1 O6 O7 size
a a lda 50 0.647 0.277 0.579 True
100 0.718 0.293 0.710 True
150 0.736 0.313 0.775 True
200 0.762 0.299 0.821 True
250 0.782 0.310 0.854 True
svm 50 0.569 0.254 0.472 True
100 0.668 0.267 0.618 True
150 0.704 0.291 0.691 True
200 0.747 0.295 0.763 True
250 0.788 0.317 0.804 True
b lda 50 0.656 0.274 0.587 True
100 0.716 0.296 0.707 True
150 0.740 0.301 0.783 True
200 0.767 0.314 0.829 True
250 0.783 0.325 0.854 True
svm 50 0.559 0.258 0.482 True
100 0.648 0.280 0.604 True
150 0.698 0.286 0.704 True
200 0.725 0.303 0.766 True
250 0.734 0.316 0.800 True
b a lda 50 0.662 0.275 0.583 True
100 0.718 0.300 0.701 True
150 0.744 0.316 0.784 True
200 0.767 0.323 0.831 True
250 0.780 0.331 0.867 True
svm 50 0.578 0.266 0.479 True
100 0.672 0.272 0.606 True
150 0.725 0.303 0.708 True
200 0.759 0.320 0.771 True
250 0.790 0.336 0.814 True
b lda 50 0.659 0.270 0.587 True
100 0.723 0.296 0.707 True
150 0.742 0.326 0.786 True
200 0.769 0.324 0.825 True
250 0.784 0.338 0.867 True
svm 50 0.576 0.253 0.481 True
100 0.664 0.272 0.612 True
150 0.702 0.308 0.707 True
200 0.727 0.322 0.762 True
250 0.739 0.336 0.809 True
J-PPC lowest at every size, per reading:
O1 O6 O7
a a lda True
svm True
b lda True
svm True
b a lda True
svm True
b lda True
svm True
Result: reproduced. At every size from 50 to 250 neurons, under every reading of O1, O6 and O7, the mean J-PPC finger-type CCGP (0.253 to 0.338) is below the means of N-PPC (0.472 to 0.867) and J-MC (0.559 to 0.790) at the same size (chance 0.2).
Blocked claims¶
fig4c-contra-discriminative, fig4e-plan-vs-move, fig5a-contra-decoding, fig5bcd-npcc-data-requirements and fig5e-planning-decoding use the alternating-cues task with a No-Go condition; fig6-online-bmi, figs6-single-array-retrospective and figs7-left-hand-bmi use the online BMI sessions. The check below lists every session's task and conditions: all 12 are the ten-finger task, with no No-Go condition and no decoder output.
print(summary[["file", "subject", "task", "trials"]].to_string())
print("conditions present in any file:", sorted(set().union(*[set(s.trials.finger) for s in S])))
print("trial table columns:", list(S[0].trials.columns))
file subject task trials 0 sub-N1_ses-2019-02-12_ecephys.nwb N1 ten-finger press task, with delay 100 1 sub-N1_ses-2019-12-20_ecephys.nwb N1 ten-finger press task, with delay 100 2 sub-P1_ses-2015-08-12_ecephys.nwb P1 ten-finger press task, with delay 100 3 sub-P1_ses-2015-08-14_ecephys.nwb P1 ten-finger press task, with delay 100 4 sub-P1_ses-2015-08-17_ecephys.nwb P1 ten-finger press task, with delay 100 5 sub-P1_ses-2015-08-21_ecephys.nwb P1 ten-finger press task, with delay 100 6 sub-P1_ses-2016-04-27_ecephys.nwb P1 ten-finger press task, with delay 100 7 sub-P1_ses-2016-05-04_ecephys.nwb P1 ten-finger press task, with delay 100 8 sub-P1_ses-2016-05-11_ecephys.nwb P1 ten-finger press task, with delay 100 9 sub-P1_ses-2016-05-16_ecephys.nwb P1 ten-finger press task, with delay 100 10 sub-P1_ses-2016-05-20_ecephys.nwb P1 ten-finger press task, with delay 100 11 sub-P1_ses-2016-05-23_ecephys.nwb P1 ten-finger press task, with delay 100 conditions present in any file: ['Li', 'Lm', 'Lp', 'Lr', 'Lt', 'Ri', 'Rm', 'Rp', 'Rr', 'Rt'] trial table columns: ['Run', 'TrialNumber', 'finger', 'hand', 'finger_type', 'start_time', 'stop_time', 'cue_on_time', 'cue_off_time', 'go_on_time']
# The dandiset this one continues (cell 3): DANDI:000147, its description, related resources and files
m147 = requests.get(f"{F.API}/dandisets/000147/versions/0.221122.2256/", timeout=60).json()
print("DANDI:000147 @ 0.221122.2256:", m147["name"])
print(m147["description"])
print("related:", [(r.get("relation"), r.get("identifier"), r.get("name")) for r in m147.get("relatedResource", [])])
a147 = requests.get(f"{F.API}/dandisets/000147/versions/0.221122.2256/assets/?page_size=100", timeout=60).json()["results"]
print(len(a147), "files:", sorted(a["path"] for a in a147))
print(f"notebook run time so far: {(time.time() - T0) / 60:.1f} min")
DANDI:000147 @ 0.221122.2256: PPC_Finger: human posterior parietal cortex recordings during attempted finger movements
This dataset contains neural activity and trial information from a tetraplegic human participant performing a brain-computer interface (BCI) finger press task. Neural activity was recorded from a 96-channel Utah array (Blackrock) implanted in the left posterior parietal cortex (PPC) at the junction of the postcentral and intraparietal sulci (PC-IP). On each trial, a finger cue, at a pseudorandom location on a screen, was selected. The participant immediately looked at the cue and pressed the cued finger until decoder feedback was shown (1.5 seconds after the cue).
This dataset includes sorted unit spiking times, finger cue, in-session classified finger, cue location, and trial timing.
related: [('dcite:IsDescribedBy', 'https://doi.org/10.7554/eLife.74478', 'Stability of motor representations after paralysis')]
10 files: ['sub-P1/sub-P1_ses-2018-09-10_ecephys.nwb', 'sub-P1/sub-P1_ses-2018-09-17a_ecephys.nwb', 'sub-P1/sub-P1_ses-2018-09-17b_ecephys.nwb', 'sub-P1/sub-P1_ses-2018-09-24_ecephys.nwb', 'sub-P1/sub-P1_ses-2018-09-26_ecephys.nwb', 'sub-P1/sub-P1_ses-2018-10-01_ecephys.nwb', 'sub-P1/sub-P1_ses-2018-10-12_ecephys.nwb', 'sub-P1/sub-P1_ses-2018-10-15_ecephys.nwb', 'sub-P1/sub-P1_ses-2018-10-17_ecephys.nwb', 'sub-P1/sub-P1_ses-2018-10-22_ecephys.nwb'] notebook run time so far: 14.3 min
Summary¶
| Claim | Outcome |
|---|---|
| tenfinger-recording-counts | reproduced (counts before the 0.5 Hz exclusion) |
| tenfinger-pooled-npcc-neurons | reproduced (1114 before exclusion) |
| tenfinger-discriminative-fraction | not reproduced (at most 48.9%, 38.1%, 65.6% against 66%, 57%, 78%; a 1 s window, not the paper's, would match) |
| tenfinger-decoding-accuracy | reproduced |
| tenfinger-error-structure | reproduced |
| fig8b-matching-closer | reproduced, resting on an assumption (noise scaling of rsatoolbox) |
| fig8a-npcc-index-distinct | reproduced |
| fig8e-hand-ccgp | reproduced |
| fig8e-finger-ccgp | reproduced (under O9 b; J-MC below the 0.9 margin of O9 a) |
| jppc-finger-ccgp-low | partially reproduced, resting on an assumption (standard accuracy on the pseudo-population; p < 0.001 in one reading of eight; per-session accuracy would give reproduced) |
| figs11-jppc-lower-matched-size | reproduced |
| fig4c, fig4e, fig5a, fig5bcd, fig5e | blocked, data_missing |
| fig6-online-bmi, figs6, figs7 | blocked, data_missing |
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 links no code. The 0.5 Hz exclusion, fold assignment, noise covariance, partitions, pseudo-trial matching and decoder of the CCGP are run under each reading the claims file names. Two choices are assumptions: the per-partition noise precision passed to rsatoolbox is that of the partition mean (the toolbox's documented input), and the standard accuracy that bounds CCGP is computed on the same pseudo-population. The CCGP null shuffles condition labels across pseudo-trials, and for random matching (O6 b) draws a new matching for each shuffle, which makes its p conservative. The standard accuracy is also shown per session under every reading.
What was found about the dataset. It holds the ten-finger task only (10 sessions of N, 2 of J), not the alternating-cues or BMI sessions the paper also reports. The dandiset it continues, DANDI:000147, holds a BCI finger-press task of participant N's PC-IP array with in-session classified finger, the task the paper attributes to ref 46 (10 files of subject P1 from 2018, described by doi 10.7554/eLife.74478; last code cell), not the alternating-cues task, and the paper does not name it. The subject identifiers are P1 for participant N and N1 for participant J, so N1 does not mean N. Most units store their spike times out of time order. Trial columns Run, TrialNumber, cue_on_time, cue_off_time and go_on_time have no description, and the Delay screen onset is only implied by cue_off_time. The per-session unit counts the paper reports match the units before the 0.5 Hz exclusion; J-MC has 139 units in its first session against the paper's maximum of 130.