Optimum-1: Human Single-Unit Auditory Deviance Detection

Analysis code and de-identified data for the study of single-unit and multi-unit responses to auditory deviants in the human brain, using the Optimum-1 multi-feature oddball paradigm.

Citation

Vinícius Carvalho, Alejandro Nasimbera, Silvia Kochen, et al. Human single-neuron responses to multi-feature auditory deviants: Evidence from amygdala, hippocampus, insula, and auditory cortex. Authorea. December 06, 2025.
DOI: 10.22541/au.176502255.57701730/v1

Repository Structure

├── MMN_MSfigs.py                    # Main analysis script (all figures and statistics)
├── plot_combined_SDF_raster_mmn.py  # Standalone SDF + raster plot generator (Figs 2, 8–9)
├── MMN_funcs.py                     # Helper functions (SDF building, GLM fitting, etc.)
├── MMN_stathelp.py                  # Nonparametric statistics (permutation tests, FDR)
├── Optimum1_ProcessedData/          # De-identified data (download from OSF; see Data section)
│   ├── Session1/
│   │   ├── SDFs/                    # P01_analysis.mat, P02_analysis.mat, ...
│   │   └── FRTabs/                  # P01_FRtable.mat, P01_FRtable.csv, ...
│   ├── Session2/
│   └── Session3/
└── Figures/                         # Precomputed model outputs (download from OSF; also created by a full run)
    └── Session123/
        ├── OutlierInfo/             # Per-subject outlier trial indices (.pkl), by session
        ├── EpochUnitAnalysis_SDF/   # Epoch GLM results
        │   └── GLM_NormRAgaussian_identity/minFR1/
        │       ├── combined/        # epoch_results.pkl, epoch_summary_table.csv, fisher_results.csv
        │       ├── category/
        │       └── individual/
        ├── LME_timeresolved_smoothminFR1/       # Time-resolved LME, all units
        │   ├── combined/            # results_{Region}_NormRA.pkl  (Amy, HPC, TTG, Ins)
        │   └── CAT+IND/
        ├── LME_timeresolved_SU_smoothminFR1/    # Time-resolved LME, single units only
        │   ├── combined/            # results_{Region}_NormRA.pkl  (Amy, HPC)
        │   └── CAT+IND/
        ├── LME_timeresolved_MU_smoothminFR1/    # Time-resolved LME, multi-units only
        │   ├── combined/
        │   └── CAT+IND/
        └── Carryover_Pairwise_LME_timeresolved_smooth_minFR1/
            ├── category/            # pairwise_results_{Region}_NormRA.pkl
            └── individual/

Usage

Option A — Regenerate figures from precomputed outputs (fast)

This skips the expensive GLM/LME fitting (which can take several hours) and loads pre-saved results directly.

1. Download both Optimum1_ProcessedData/ and Figures/ from OSF (https://osf.io/exj6r/) and place them next to the scripts.
2. Open MMN_MSfigs.py and confirm the load flags are set to True (they are by default):

   loadResultsEpoch = True   # ~line 2134 — skips epoch GLM fitting
   # loadResults=True is the default in each run_mmn_glm_analysis / run_carryover_pairwise_analysis call

3. Run cells sequentially in Spyder or VS Code. The data loading section (building SDF DataFrames from the .mat files) still runs as normal; only the model-fitting steps are skipped.

Option B — Full reanalysis from scratch

1. Download Optimum1_ProcessedData/ from OSF and set the data path:

   dataPath = Path('./Optimum1_ProcessedData')

2. Set all load flags to False:

   loadResultsEpoch = False   # ~line 2134
   # set loadResults=False in each run_mmn_glm_analysis / run_carryover_pairwise_analysis call

3. Run cells sequentially. The time-resolved LME sections are parallelised with joblib but remain computationally intensive.

Running the Raster Plot Script

# Edit configuration in plot_combined_SDF_raster_mmn.py, then:
python plot_combined_SDF_raster_mmn.py

Manuscript Figures

Figure	Description	Script / Section
Fig 1	Task design and electrode recording locations	schematic + MNI panel; not code-generated
Fig 2	Representative SDF + raster plots — Amy and HPC units	plot_combined_SDF_raster_mmn.py
Fig 3	Dataset overview: unit counts, firing rates, participant/session distribution	MMN_MSfigs.py — Main execution
Figs 4–6	Epoch-based LME/GLM — combined, category, and individual deviant models; heatmap + raincloud panels (Amy, HPC)	MMN_MSfigs.py — #%% Build 3 manuscript epoch figures (heatmap + raincloud combined)
Fig 4b	Per-unit response heatmap by participant — category + combined model	MMN_MSfigs.py — #%% Per-Unit Heatmap Figure — Category + Combined, REFINED
Figs 7, S2–S4	Time-resolved LME — combined and CAT+IND models, Amy/HPC and TTG/PIns	MMN_MSfigs.py — #%% MS region-pair figures (time-resolved ensemble analysis)
Figs 8–9	TTG and PIns unit responses — SDF + raster plots	plot_combined_SDF_raster_mmn.py
Fig S1	T1 MRI slices with estimated electrode tip locations (TTG and PIns probes)	not code-generated
Fig S5–S6	Carryover effects — time-resolved individual deviants (Amy/HPC and TTG/PIns)	MMN_MSfigs.py — #%% MS region-pair figures (carryover analysis)
Fig S7	Time-resolved LME — SU and MU separately (Amy, HPC)	MMN_MSfigs.py — #%% MS region-pair figures (time-resolved ensemble analysis), SU/MU
Supp Table S2	Per-patient summary: unit counts, SU/MU, firing rates, modulation rates	MMN_MSfigs.py — #%% Supp Table: per-patient × brain-region summary

Analysis Overview

Statistical Approach

Unit-level epoch analyses use GLMs. Population-level (ensemble) time-resolved analyses use linear mixed-effects models (LME) with participant as a random intercept (model_type='LME'). GLM infrastructure is retained in the code for comparison.

Analysis Pipeline

1. Data loading - Load de-identified .mat files, build SDF DataFrames, apply rolling-average normalization
2. Unit overview (Fig 3) - Firing rates, unit counts, participant distribution
3. Epoch-based analysis (Figs 4–6) - GLM per unit within Early/Middle/Late epochs; population-level Fisher combined test with FDR correction
4. Time-resolved analysis (Figs 7, S2–S4) - LME at 1ms resolution with delta-AIC model comparison; FDR-corrected significance
5. Carryover analysis (Figs S5–S6) - Pairwise analysis of standard responses conditioned on preceding deviant type
6. SDF + raster plots (Figs 2, 8–9) - Individual unit plots with sliding-window permutation statistics

Key Parameters

• Normalization: Rolling-average (sdf_norm='RA'), period: -100 to 450 ms
• Decimation: 10x (from 1kHz to 100Hz time resolution)
• Minimum firing rate: 1 Hz
• Outlier removal: 4x IQR multiplier
• Permutation tests: 2000-5000 permutations, FDR alpha = 0.05

Data

De-identified data are available on OSF: https://osf.io/exj6r/

De-identified Data Format

Each participant file (P01_analysis.mat) contains:

• sdf - Spike density functions (trial x unit x time)
• spikeTrials - Raw spike times per trial
• trlInfo - Trial metadata (event IDs, block numbers)
• ALLspikesInfo - Unit metadata (brain region, firing rate, SUA/MUA classification)

13 participants (P01-P13), 17 sessions total. Three participants have multiple sessions:

• P02: Sessions 1, 2
• P06: Sessions 1, 3
• P08: Sessions 1, 2

Brain Regions

• HPC - Hippocampus (46 units)
• Amy - Amygdala (48 units)
• TTG - Transverse Temporal Gyrus (2 units)
• PIns - Posterior Insula (2 units)

Requirements

numpy
pandas
matplotlib
seaborn
scipy
statsmodels
pymatreader
joblib