0
# Statistical Analysis
1

2
Statistical methods including cluster-based permutation testing, multiple comparisons correction, and specialized analysis tools designed for neuroimaging data with high temporal and spatial resolution.
3

4
## Capabilities
5

6
### Multiple Comparisons Correction
7

8
Correct for multiple comparisons across time, frequency, and spatial dimensions.
9

10
```python { .api }
11
def bonferroni_correction(pvals: ArrayLike, alpha: float = 0.05) -> Tuple[ArrayLike, float]:
12
    """
13
    Bonferroni correction for multiple comparisons.
14
    
15
    Parameters:
16
    - pvals: P-values to correct
17
    - alpha: Family-wise error rate
18
    
19
    Returns:
20
    Tuple of (reject array, corrected alpha level)
21
    """
22

23
def fdr_correction(pvals: ArrayLike, alpha: float = 0.05, method: str = 'indep') -> Tuple[ArrayLike, ArrayLike]:
24
    """
25
    False discovery rate correction.
26
    
27
    Parameters:
28
    - pvals: P-values to correct
29
    - alpha: False discovery rate
30
    - method: FDR method ('indep' or 'negcorr')
31
    
32
    Returns:
33
    Tuple of (reject array, corrected p-values)
34
    """
35
```
36

37
### Cluster-Based Permutation Testing
38

39
Non-parametric statistical testing that accounts for multiple comparisons by clustering adjacent significant tests.
40

41
```python { .api }
42
def permutation_cluster_test(X: ArrayLike, threshold: Optional[float] = None, n_permutations: int = 1024,
43
                            tail: int = 0, stat_fun: Optional[callable] = None,
44
                            adjacency: Optional[ArrayLike] = None, n_jobs: int = 1,
45
                            seed: Optional[int] = None, max_step: int = 1, exclude: Optional[ArrayLike] = None,
46
                            step_down_p: float = 0, t_power: float = 1, out_type: str = 'indices',
47
                            check_disjoint: bool = False, buffer_size: Optional[int] = None,
48
                            verbose: Optional[Union[bool, str, int]] = None) -> Tuple:
49
    """
50
    Cluster-based permutation test for comparing conditions.
51
    
52
    Parameters:
53
    - X: Data arrays to compare (list of arrays or single array)
54
    - threshold: Statistical threshold for clustering
55
    - n_permutations: Number of permutations
56
    - tail: Test direction (0=two-tailed, 1=greater, -1=less)
57
    - stat_fun: Statistical function to use
58
    - adjacency: Adjacency matrix for clustering
59
    - n_jobs: Number of parallel jobs
60
    - seed: Random seed for reproducibility
61
    - max_step: Maximum cluster step size
62
    - exclude: Samples to exclude from clustering
63
    - step_down_p: Step-down p-value threshold
64
    - t_power: Power for threshold-free cluster enhancement
65
    - out_type: Output type ('indices' or 'mask')
66
    - check_disjoint: Check for disjoint clusters
67
    - buffer_size: Buffer size for permutations
68
    - verbose: Verbosity level
69
    
70
    Returns:
71
    Tuple of (T_obs, clusters, cluster_p_values, H0)
72
    """
73

74
def permutation_cluster_1samp_test(X: ArrayLike, threshold: Optional[float] = None, n_permutations: int = 1024,
75
                                  tail: int = 0, stat_fun: Optional[callable] = None,
76
                                  adjacency: Optional[ArrayLike] = None, n_jobs: int = 1,
77
                                  seed: Optional[int] = None, max_step: int = 1,
78
                                  exclude: Optional[ArrayLike] = None, step_down_p: float = 0,
79
                                  t_power: float = 1, out_type: str = 'indices',
80
                                  check_disjoint: bool = False, buffer_size: Optional[int] = None,
81
                                  verbose: Optional[Union[bool, str, int]] = None) -> Tuple:
82
    """
83
    One-sample cluster-based permutation test.
84
    
85
    Parameters:
86
    - X: Data array (n_observations, n_features)
87
    - threshold: Statistical threshold for clustering
88
    - n_permutations: Number of permutations
89
    - tail: Test direction
90
    - stat_fun: Statistical function
91
    - adjacency: Adjacency matrix for clustering
92
    - n_jobs: Number of parallel jobs
93
    - seed: Random seed
94
    - max_step: Maximum cluster step size
95
    - exclude: Samples to exclude
96
    - step_down_p: Step-down p-value threshold
97
    - t_power: Power for TFCE
98
    - out_type: Output type
99
    - check_disjoint: Check for disjoint clusters
100
    - buffer_size: Buffer size
101
    - verbose: Verbosity level
102
    
103
    Returns:
104
    Tuple of (T_obs, clusters, cluster_p_values, H0)
105
    """
106

107
def spatio_temporal_cluster_test(X: ArrayLike, adjacency: ArrayLike, threshold: Optional[float] = None,
108
                                n_permutations: int = 1024, tail: int = 0, stat_fun: Optional[callable] = None,
109
                                n_jobs: int = 1, seed: Optional[int] = None, max_step: int = 1,
110
                                exclude: Optional[ArrayLike] = None, step_down_p: float = 0,
111
                                t_power: float = 1, out_type: str = 'indices',
112
                                check_disjoint: bool = False, buffer_size: Optional[int] = None,
113
                                verbose: Optional[Union[bool, str, int]] = None) -> Tuple:
114
    """
115
    Spatio-temporal cluster-based permutation test.
116
    
117
    Parameters:
118
    - X: Data arrays (list of n_conditions arrays of shape n_observations × n_vertices × n_times)
119
    - adjacency: Spatial adjacency matrix
120
    - threshold: Statistical threshold
121
    - n_permutations: Number of permutations
122
    - tail: Test direction
123
    - stat_fun: Statistical function
124
    - n_jobs: Number of parallel jobs
125
    - seed: Random seed
126
    - max_step: Maximum cluster step size
127
    - exclude: Samples to exclude
128
    - step_down_p: Step-down threshold
129
    - t_power: TFCE power
130
    - out_type: Output type
131
    - check_disjoint: Check disjointness
132
    - buffer_size: Buffer size
133
    - verbose: Verbosity level
134
    
135
    Returns:
136
    Tuple of (T_obs, clusters, cluster_p_values, H0)
137
    """
138

139
def spatio_temporal_cluster_1samp_test(X: ArrayLike, adjacency: ArrayLike, threshold: Optional[float] = None,
140
                                      n_permutations: int = 1024, tail: int = 0,
141
                                      stat_fun: Optional[callable] = None, n_jobs: int = 1,
142
                                      seed: Optional[int] = None, max_step: int = 1,
143
                                      exclude: Optional[ArrayLike] = None, step_down_p: float = 0,
144
                                      t_power: float = 1, out_type: str = 'indices',
145
                                      check_disjoint: bool = False, buffer_size: Optional[int] = None,
146
                                      verbose: Optional[Union[bool, str, int]] = None) -> Tuple:
147
    """
148
    One-sample spatio-temporal cluster test.
149
    
150
    Parameters:
151
    - X: Data array (n_observations × n_vertices × n_times)
152
    - adjacency: Spatial adjacency matrix
153
    - threshold: Statistical threshold
154
    - n_permutations: Number of permutations
155
    - tail: Test direction
156
    - stat_fun: Statistical function
157
    - n_jobs: Number of parallel jobs
158
    - seed: Random seed
159
    - max_step: Maximum cluster step size
160
    - exclude: Samples to exclude
161
    - step_down_p: Step-down threshold
162
    - t_power: TFCE power
163
    - out_type: Output type
164
    - check_disjoint: Check disjointness
165
    - buffer_size: Buffer size
166
    - verbose: Verbosity level
167
    
168
    Returns:
169
    Tuple of (T_obs, clusters, cluster_p_values, H0)
170
    """
171
```
172

173
### Parametric Statistical Tests
174

175
Classical parametric tests adapted for neuroimaging data structures.
176

177
```python { .api }
178
def f_oneway(*args: ArrayLike) -> Tuple[ArrayLike, ArrayLike]:
179
    """
180
    One-way ANOVA F-test.
181
    
182
    Parameters:
183
    - *args: Data arrays for each group
184
    
185
    Returns:
186
    Tuple of (F-statistics, p-values)
187
    """
188

189
def f_mway_rm(data: ArrayLike, factor_levels: ArrayLike, effects: str = 'A*B',
190
             correction: bool = False, alpha: float = 0.05) -> Dict:
191
    """
192
    Multi-way repeated measures ANOVA.
193
    
194
    Parameters:
195
    - data: Data array (n_subjects × n_conditions × n_times)
196
    - factor_levels: Number of levels for each factor
197
    - effects: Effects to test ('A', 'B', 'A*B', etc.)
198
    - correction: Apply sphericity correction
199
    - alpha: Significance level
200
    
201
    Returns:
202
    Dictionary with F-statistics, p-values, and effect information
203
    """
204

205
def ttest_1samp_no_p(X: ArrayLike, popmean: float = 0, sigma: Optional[float] = None) -> ArrayLike:
206
    """
207
    One-sample t-test without p-value computation (for speed).
208
    
209
    Parameters:
210
    - X: Data array
211
    - popmean: Population mean to test against
212
    - sigma: Standard deviation (if known)
213
    
214
    Returns:
215
    T-statistics array
216
    """
217

218
def ttest_ind_no_p(X1: ArrayLike, X2: ArrayLike, equal_var: bool = True) -> ArrayLike:
219
    """
220
    Independent samples t-test without p-values.
221
    
222
    Parameters:
223
    - X1: First group data
224
    - X2: Second group data  
225
    - equal_var: Assume equal variances
226
    
227
    Returns:
228
    T-statistics array
229
    """
230
```
231

232
### Non-Parametric Tests and Bootstrap
233

234
Non-parametric alternatives and bootstrap methods for statistical inference.
235

236
```python { .api }
237
def permutation_t_test(X: ArrayLike, n_permutations: int = 10000, tail: int = 0,
238
                      n_jobs: int = 1, seed: Optional[int] = None,
239
                      verbose: Optional[Union[bool, str, int]] = None) -> Tuple[ArrayLike, ArrayLike]:
240
    """
241
    Permutation-based t-test.
242
    
243
    Parameters:
244
    - X: Data array (n_observations, n_features)
245
    - n_permutations: Number of permutations
246
    - tail: Test direction (0=two-tailed, 1=greater, -1=less)
247
    - n_jobs: Number of parallel jobs
248
    - seed: Random seed
249
    - verbose: Verbosity level
250
    
251
    Returns:
252
    Tuple of (T-statistics, p-values)
253
    """
254

255
def bootstrap_confidence_interval(arr: ArrayLike, ci: float = 0.95, n_bootstraps: int = 10000,
256
                                 stat_fun: str = 'mean', seed: Optional[int] = None,
257
                                 n_jobs: int = 1, verbose: Optional[Union[bool, str, int]] = None) -> ArrayLike:
258
    """
259
    Bootstrap confidence intervals.
260
    
261
    Parameters:
262
    - arr: Data array
263
    - ci: Confidence interval level
264
    - n_bootstraps: Number of bootstrap samples
265
    - stat_fun: Statistic to compute ('mean', 'median', 'std')
266
    - seed: Random seed
267
    - n_jobs: Number of parallel jobs
268
    - verbose: Verbosity level
269
    
270
    Returns:
271
    Confidence interval bounds
272
    """
273
```
274

275
### Linear Regression Analysis
276

277
Regression methods for relating neural signals to experimental variables.
278

279
```python { .api }
280
def linear_regression(design_matrix: ArrayLike, data: ArrayLike, names: Optional[List[str]] = None) -> Dict:
281
    """
282
    Compute linear regression.
283
    
284
    Parameters:
285
    - design_matrix: Design matrix (n_observations × n_regressors)
286
    - data: Data matrix (n_observations × n_features)
287
    - names: Regressor names
288
    
289
    Returns:
290
    Dictionary with regression coefficients, t-statistics, p-values, etc.
291
    """
292

293
def linear_regression_raw(raw: Raw, design_matrix: ArrayLike, names: Optional[List[str]] = None,
294
                         tmin: float = 0.0, tmax: Optional[float] = None,
295
                         reject_by_annotation: bool = True, picks: Optional[Union[str, List]] = None) -> Dict:
296
    """
297
    Linear regression on raw data.
298
    
299
    Parameters:
300
    - raw: Raw data instance
301
    - design_matrix: Design matrix
302
    - names: Regressor names
303
    - tmin: Start time
304
    - tmax: End time
305
    - reject_by_annotation: Reject annotated segments
306
    - picks: Channel selection
307
    
308
    Returns:
309
    Dictionary with regression results
310
    """
311
```
312

313
### Adjacency Matrices and Clustering Utilities
314

315
Construct spatial adjacency matrices for cluster-based statistics.
316

317
```python { .api }
318
def combine_adjacency(adjacency_list: List[ArrayLike]) -> ArrayLike:
319
    """
320
    Combine multiple adjacency matrices.
321
    
322
    Parameters:
323
    - adjacency_list: List of adjacency matrices
324
    
325
    Returns:
326
    Combined adjacency matrix
327
    """
328

329
def spatial_src_adjacency(src: SourceSpaces, verbose: Optional[Union[bool, str, int]] = None) -> ArrayLike:
330
    """
331
    Compute adjacency matrix for source space.
332
    
333
    Parameters:
334
    - src: Source space
335
    - verbose: Verbosity level
336
    
337
    Returns:
338
    Source space adjacency matrix
339
    """
340

341
def spatial_tris_adjacency(tris: ArrayLike, remap_vertices: bool = False) -> ArrayLike:
342
    """
343
    Compute adjacency matrix from triangulation.
344
    
345
    Parameters:
346
    - tris: Triangle array (n_tris × 3)
347
    - remap_vertices: Remap vertex indices
348
    
349
    Returns:
350
    Spatial adjacency matrix
351
    """
352

353
def spatial_dist_adjacency(src: SourceSpaces, dist: float, verbose: Optional[Union[bool, str, int]] = None) -> ArrayLike:
354
    """
355
    Compute distance-based adjacency matrix.
356
    
357
    Parameters:
358
    - src: Source space with distance information
359
    - dist: Maximum distance for adjacency
360
    - verbose: Verbosity level
361
    
362
    Returns:
363
    Distance-based adjacency matrix
364
    """
365

366
def spatial_inter_hemi_adjacency(src: SourceSpaces, dist: float,
367
                                verbose: Optional[Union[bool, str, int]] = None) -> ArrayLike:
368
    """
369
    Get inter-hemisphere adjacency.
370
    
371
    Parameters:
372
    - src: Source space
373
    - dist: Maximum inter-hemisphere distance
374
    - verbose: Verbosity level
375
    
376
    Returns:
377
    Inter-hemisphere adjacency matrix
378
    """
379
```
380

381
### Cluster Analysis and Summarization
382

383
Analyze and summarize results from cluster-based permutation tests.
384

385
```python { .api }
386
def summarize_clusters_stc(clu: Tuple, tstep: float, tmin: float, subject: str,
387
                          vertices: Optional[List[ArrayLike]] = None,
388
                          verbose: Optional[Union[bool, str, int]] = None) -> List[Dict]:
389
    """
390
    Summarize cluster test results in source space.
391
    
392
    Parameters:
393
    - clu: Cluster test results tuple
394
    - tstep: Time step
395
    - tmin: Start time
396
    - subject: Subject name
397
    - vertices: Vertex arrays (if different from clu)
398
    - verbose: Verbosity level
399
    
400
    Returns:
401
    List of cluster summary dictionaries
402
    """
403
```
404

405
## Usage Examples
406

407
### Basic Cluster-Based Permutation Test
408

409
```python
410
import mne
411
from mne.stats import permutation_cluster_test
412
import numpy as np
413

414
# Load epochs for two conditions
415
epochs1 = mne.read_epochs('condition1-epo.fiv')  
416
epochs2 = mne.read_epochs('condition2-epo.fiv')
417

418
# Get data arrays
419
X = [epochs1.get_data(), epochs2.get_data()]
420

421
# Run cluster-based permutation test
422
T_obs, clusters, cluster_p_values, H0 = permutation_cluster_test(
423
    X, n_permutations=1000, threshold=None, tail=0, n_jobs=1)
424

425
# Find significant clusters
426
significant_clusters = [c for c, p in zip(clusters, cluster_p_values) if p < 0.05]
427
print(f"Found {len(significant_clusters)} significant clusters")
428
```
429

430
### Spatio-Temporal Cluster Test
431

432
```python
433
import mne
434
from mne.stats import spatio_temporal_cluster_1samp_test
435

436
# Load epochs
437
epochs = mne.read_epochs('sample-epo.fiv')
438

439
# Get data (subtract baseline)
440
X = epochs.get_data() - epochs.get_data().mean(axis=-1, keepdims=True)
441

442
# Create adjacency matrix for channels
443
adjacency, ch_names = mne.channels.find_ch_adjacency(epochs.info, ch_type='eeg')
444

445
# Run spatio-temporal cluster test  
446
T_obs, clusters, cluster_p_values, H0 = spatio_temporal_cluster_1samp_test(
447
    X, adjacency=adjacency, n_permutations=1000, threshold=None, tail=1)
448

449
# Visualize results
450
evoked = mne.EvokedArray(T_obs, epochs.info, tmin=epochs.tmin)
451
evoked.plot_topomap(times=np.linspace(0.05, 0.3, 6))
452
```
453

454
### Source Space Cluster Test
455

456
```python
457
import mne
458
from mne.stats import spatio_temporal_cluster_1samp_test, summarize_clusters_stc
459

460
# Load source estimates for multiple subjects
461
stcs = []
462
for subject in subjects:
463
    stc = mne.read_source_estimate(f'{subject}-stc.fiv')
464
    stcs.append(stc)
465

466
# Convert to array (n_subjects × n_vertices × n_times)
467
X = np.array([stc.data for stc in stcs])
468

469
# Get source space adjacency
470
src = mne.read_source_spaces('fsaverage-oct6-src.fiv')
471
adjacency = mne.spatial_src_adjacency(src)
472

473
# Run cluster test
474
T_obs, clusters, cluster_p_values, H0 = spatio_temporal_cluster_1samp_test(
475
    X, adjacency=adjacency, n_permutations=1000, threshold=None, tail=1)
476

477
# Summarize significant clusters
478
significant_clusters = summarize_clusters_stc(
479
    (T_obs, clusters, cluster_p_values, H0), 
480
    tstep=stcs[0].tstep, tmin=stcs[0].tmin, subject='fsaverage')
481

482
# Create source estimate with T-statistics
483
stc_T = mne.SourceEstimate(T_obs, vertices=stcs[0].vertices, 
484
                          tmin=stcs[0].tmin, tstep=stcs[0].tstep,
485
                          subject='fsaverage')
486
stc_T.plot(subject='fsaverage', subjects_dir=subjects_dir)
487
```
488

489
### Multiple Comparisons Correction
490

491
```python
492
import mne
493
from mne.stats import fdr_correction, bonferroni_correction
494
import numpy as np
495

496
# Simulate p-values from multiple tests
497
n_tests = 1000
498
p_values = np.random.uniform(0, 1, n_tests)
499
p_values[:50] = np.random.uniform(0, 0.01, 50)  # Add some "significant" values
500

501
# Apply FDR correction
502
reject_fdr, p_fdr = fdr_correction(p_values, alpha=0.05)
503
print(f"FDR: {np.sum(reject_fdr)} significant tests")
504

505
# Apply Bonferroni correction  
506
reject_bonf, alpha_bonf = bonferroni_correction(p_values, alpha=0.05)
507
print(f"Bonferroni: {np.sum(reject_bonf)} significant tests")
508
print(f"Bonferroni corrected alpha: {alpha_bonf}")
509
```
510

511
### Linear Regression Analysis
512

513
```python
514
import mne
515
from mne.stats import linear_regression
516
import numpy as np
517

518
# Load epochs
519
epochs = mne.read_epochs('sample-epo.fiv')
520

521
# Create design matrix (e.g., stimulus intensity)
522
n_epochs = len(epochs)
523
stimulus_intensity = np.random.randn(n_epochs)  
524
design_matrix = np.column_stack([np.ones(n_epochs), stimulus_intensity])
525

526
# Get data
527
data = epochs.get_data().reshape(n_epochs, -1)  # Flatten channels and time
528

529
# Run regression
530
results = linear_regression(design_matrix, data, 
531
                          names=['Intercept', 'Stimulus_Intensity'])
532

533
# Extract results
534
beta = results['beta']  # Regression coefficients
535
t_val = results['t_val']  # T-statistics  
536
p_val = results['p_val']  # P-values
537

538
# Reshape back to original dimensions
539
n_channels, n_times = epochs.get_data().shape[1:]
540
t_val_reshaped = t_val[1].reshape(n_channels, n_times)  # Stimulus effect
541

542
# Create evoked with T-statistics
543
evoked_t = mne.EvokedArray(t_val_reshaped, epochs.info, tmin=epochs.tmin)
544
evoked_t.plot_topomap(times=[0.1, 0.2, 0.3])
545
```
546

547
## Types
548

549
```python { .api }
550
import numpy as np
551
from typing import Union, Optional, List, Dict, Tuple, Any, Callable
552

553
ArrayLike = Union[np.ndarray, List, Tuple]
554
```